ebcot_encoder.c

关于视频压缩的jpeg2000压缩算法,C编写

/*****************************************************************************/
/* Copyright 1998, Hewlett-Packard Company                                   */
/* All rights reserved                                                       */
/* File: "ebcot_encoder.h"                                                   */
/* Description: Central source file for EBCOT encoder                        */
/* Author: David Taubman                                                     */
/* Affiliation: Hewlett-Packard and                                          */
/*              The University of New South Wales, Australia                 */
/* Version: VM6.0                                                            */
/* Last Revised: 28 January, 2000                                            */
/*****************************************************************************/

/*****************************************************************************/
/* Modified to incorporate MQ-coder by Mitsubishi Electric Corp.             */
/* Copyright 1999, Mitsubishi Electric Corp.                                 */
/* All rights reserved for modified parts                                    */
/*****************************************************************************/

/*****************************************************************************/
/* Modified to incorporate optional entropy coder                            */
/* Copyright 1999 Science Applications International Corporation (SAIC).     */
/* All Rights Reservedi for modified parts.                                  */
/*****************************************************************************/

/*****************************************************************************/
/* Modifications tagged with comment lines or delimiters involving the       */
/* string, "David T mod" or "mod by David T", have been made by David        */
/* Taubman; they are copyrighted by Hewlett-Packard Company with all rights  */
/* reserved for the modified parts.                                          */
/*****************************************************************************/

/*****************************************************************************/
/* Modified by David Taubman to support new termination methods for parallel */
/* coding options, to enable more efficient implementations and to support   */
/* error resilience.  Copyright 1999 by Hewlett-Packard Company with all     */
/* rights reserved for the modified parts.                                   */
/*****************************************************************************/

/*****************************************************************************/
/* Modified to combine entropy coders                                        */
/* Copyright 1999 Science Applications International Corporation (SAIC).     */
/* Copyright 1999 University of Arizona, Arizona Board of Regents.           */
/* All Rights Reservedi for modified parts.                                  */
/*****************************************************************************/

/*****************************************************************************/
/* Modified to include TCQ                                                   */
/* Copyright 1999 Science Applications International Corporation (SAIC).     */
/* Copyright 1995 University of Arizona, Arizona Board of Regents.           */
/* All Rights Reservedi for modified parts.                                  */
/*****************************************************************************/

/*****************************************************************************/
/* Modified by David Taubman to include masking-sensitive distortion         */
/* calculations for R-D optimization (`-Cvis' option).  Material identified  */
/* by "David T Cvis mod" comments has been added by David Taubman; it is     */
/* copyrighted by the University of New South Wales (Copyright 1999) with    */
/* all rights reserved for the modified parts.                               */
/*****************************************************************************/

/*****************************************************************************/
/* Modified by David Taubman to support arbitrary reference points for the   */
/* transform and the various regular partitions, so as to facilitate         */
/* cropping and geometric transformations in the compressed domain and to    */
/* enable full support for CRF's single-sample overlap Wavelet transform,    */
/* as originally documented in Seoul.  Changes are too numerous to flag      */
/* individually within the code.  Changes copyrighted by HP with all rights  */
/* reserved for the modified parts.                                          */
/*****************************************************************************/

/*****************************************************************************/
/* Modified by David Taubman to support interface modifications, arbitrary   */
/* changes in coding parameters from component to component and from tile    */
/* to tile, packet partitions, rich packet sequencing conventions and to     */
/* support the full generality of PART-1 of the JPEG2000                     */
/* standard, and to support most anticipated generality of PART-2.  Changes  */
/* are too numerous to flag individually within the code, which in some      */
/* places has been completely rewritten.  All changes copyrighted by HP with */
/* all rights reserved for the modified parts.                               */
/*****************************************************************************/

#include <local_services.h>
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <math.h>
#include <time.h>
#include <assert.h>
#include <ifc.h>
#include <markers.h>
#include "ebcot_common.h"
#include "ebcot_encoder.h"
#include "ebcot_send_bits.h"
#include "ebcot_vpw.h"


/* ========================================================================= */
/* --------------------------- Utility Functions --------------------------- */
/* ========================================================================= */

/*****************************************************************************/
/* STATIC                       is_power_of_2                                */
/*****************************************************************************/

static int
  is_power_of_2(int val)
{
  if (val == 0)
    return(0);
  while (!(val & 1))
    val>>=1;
  return(val == 1);
}

/*****************************************************************************/
/* STATIC                      start_symbol_log                              */
/*****************************************************************************/
#ifdef DST_LOG_SYMBOLS
static FILE *
  start_symbol_log(ebcot_band_info_ptr band, ebcot_block_info_ptr block)
{
  FILE *fp;
  char fname[81];

  if (band->tnum != 0)
    local_error("Symbol logging is permitted only on single tile images!");
  sprintf(fname,"s_%d_%d_%d",band->level_idx,band->band_idx,
          (int)(block-band->blocks));
  fp = fopen(fname,"wb");
  assert(fp != NULL);
  return(fp);
}
#endif /* DST_LOG_SYMBOLS */

/*****************************************************************************/
/* STATIC              convert_dimension_from_level_to_band                  */
/*****************************************************************************/

static int
  convert_dimension_from_level_to_band(int dim, int is_vertical,
                                       int hp_descent, int band_idx)
{
  int d, orient, v;

  v = (is_vertical)?1:0;
  if (hp_descent == 0)
    assert(band_idx == LL_BAND);
  for (d=hp_descent-1; d >= 0; d--)
    {
      assert(band_idx != LL_BAND);
      orient = (band_idx >> (d+d)) & 3;
      switch (orient) {
        case LL_BAND: dim = (dim+1)>>1; break;
        case HL_BAND: dim = (dim+v)>>1; break;
        case LH_BAND: dim = (dim+1-v)>>1; break;
        case HH_BAND: dim = dim >> 1; break;
        default: assert(0);
      }
    }
  return(dim);
}

/*****************************************************************************/
/* STATIC                     allocate_pass_info                             */
/*****************************************************************************/

static void
  allocate_pass_info(ebcot_band_info_ptr band)

 /* This function allocates the `pass_info' array for the subband and
    sets the `passes' field of each code block in the subband to point to
    the relevant portion of this array.  To do this, the function invokes
    the `forward_roi__get_max_boost' function, as necessary, to determine
    the maximum number of passes which might have to be generated for
    any given code block. */

{
  ebcot_pass_info_ptr passes;
  ebcot_block_info_ptr block;
  int block_aligned, max_passes, max_planes, block_idx;

  assert(band->roi != NULL);
  band->max_roi_boost =
    band->roi->get_max_boost(band->roi,band->component_idx,band->level_idx,
                             band->band_idx,&block_aligned,
                             &(band->roi_implicit));
  if (band->max_roi_boost == 0)
    band->roi = NULL; /* No need to call `forward_roi__check_roi'. */
  if (block_aligned)
    band->max_roi_boost = 0;

  max_planes = IMPLEMENTATION_PRECISION - 1 - band->extra_lsbs;
  max_planes += band->max_roi_boost;
  if (max_planes > (IMPLEMENTATION_PRECISION-2))
    max_planes = IMPLEMENTATION_PRECISION-2;
  max_passes = 1 + PASSES_PER_BITPLANE*(max_planes-1);

  /* Sharp Labs change below */
  if ( max_planes < 0 )
    max_planes = 0;
  if ( max_passes < 0 )
    max_passes = 0;

  band->max_passes = max_passes;
  passes = (ebcot_pass_info_ptr)
    local_malloc(EBCOT_MEM_KEY,
                 sizeof(ebcot_pass_info)*(band->total_blocks*max_passes));
  for (block_idx=0; block_idx < band->total_blocks; block_idx++)
    {
      block = band->blocks + block_idx;
      block->max_passes = 0; /* Fill this in at encode time. */
      block->passes = passes;
      passes += max_passes;
    }
}


/* ========================================================================= */
/* ---------------------- Rate-Distortion Functions ------------------------ */
/* ========================================================================= */

/*****************************************************************************/
/* STATIC               convert_double_to_rd_slope_type                      */
/*****************************************************************************/

static rd_slope_type
  convert_double_to_rd_slope_type(double slope)

 /* Converts a real-valued rate-distortion slope to the exponent-mantissa
    representation associated with the `rd_slope_type' data-type, as
    explained with its definition in "ebcot_encoder.h".  Excessively
    small or large slope values are truncated to the minimum and maximum
    legal slopes, respectively. */

{
  std_short exponent, max_exponent, mantissa, max_mantissa;
  double boundary;

  max_exponent = (1<<RD_SLOPE_EXPONENT_BITS) - 1;
  max_mantissa = (1<<RD_SLOPE_MANTISSA_BITS) - 1;
  boundary = 2.0 - 0.5 / ((double)(max_mantissa+1));
    /* Largest value which should be represented by an exponent of 0. */
  for (exponent=0; exponent < (max_exponent-1); exponent++, slope*=0.5)
    if (slope <= boundary)
      break;
  if (slope > boundary)
    slope = boundary;
  mantissa = (std_short) floor(0.5 + (slope-1.0)*(double)(max_mantissa+1));
  if (mantissa < 1)
    mantissa = 1;
  if (mantissa > max_mantissa)
    mantissa = max_mantissa;
  return((exponent<<RD_SLOPE_MANTISSA_BITS) + mantissa);
}

/*****************************************************************************/
/* STATIC                     compute_rd_slopes                              */
/*****************************************************************************/

static void
  compute_rd_slopes(int num_points, double *cumulative_wmse_reduction,
                    ebcot_pass_info_ptr passes, std_int *rd_slope_rates)
{
  int n, last_n, completed_points;
  double delta_wmse, delta_bytes, last_slope, slope;

  completed_points = 0;
  last_slope = 0; /* Prevent warnings from compiler. */
  while (completed_points < num_points)
    {
      for (last_n=-1, n=0; n < num_points; n++)
        {
          delta_wmse = cumulative_wmse_reduction[n];
          delta_bytes = (double)(passes[n].cumulative_bytes);
          if (last_n >= 0)
            {
              delta_wmse -= cumulative_wmse_reduction[last_n];
              delta_bytes -= (double)(passes[last_n].cumulative_bytes);
            }
          if (delta_wmse <= 0.0)
            {
              passes[n].rd_slope = -1; /* Never use this pass. */
              if (n >= completed_points)
                completed_points = n+1;
              continue;
            }
          if ((n < completed_points) && (passes[n].rd_slope <= 0))
            continue; /* Pass does not lie on convex hull. */
          if (delta_bytes == 0.0)
            { /* Slope is infinite.  Last point cannot lie on convex hull. */
              assert(last_n >= 0);
              passes[last_n].rd_slope = 0;
              break; /* Start over again. */
            }
          slope = delta_wmse / delta_bytes;
          if ((last_n >= 0) && (slope >= last_slope))
            { /* Last point cannot lie on convex hull. */
              passes[last_n].rd_slope = 0;
              break; /* Start over again. */
            }
          last_n = n;
          last_slope = slope;
          slope *= (double)(1<<30);
          passes[n].rd_slope = convert_double_to_rd_slope_type(slope);
          if (n >= completed_points)
            completed_points = n+1;
        }
    }

  { /* Accumulate slope rates information. */
    std_int num_bytes;
    rd_slope_type slope_val, exponent;

    num_bytes = 0;
    for (n=0, exponent=(1<<RD_SLOPE_EXPONENT_BITS)-1;
         exponent >= 0; exponent--)
      {
        for (; n < num_points; n++)
          {
            slope_val = passes[n].rd_slope;
            if ((slope_val >> RD_SLOPE_MANTISSA_BITS) >= exponent)
              num_bytes = passes[n].cumulative_bytes;
            else if (slope_val > 0)
              break;
          }
        rd_slope_rates[exponent] += num_bytes;
      }
  }
}

/* ========================================================================= */
/* --------------------------- Encoding Functions -------------------------- */
/* ========================================================================= */

/*****************************************************************************/