ebcot_encoder.c

JPEG2000 EBCOT算法源码

      double log_r1, log_r2;
      int diff, m;

      for (m=n+1; !self->layer_info[m].optimize; m++);
      diff = m-n;
      if (diff == 1)
        continue;
      log_r1 = log((double) self->layer_info[n].max_cumulative_bytes);
      log_r2 = log((double) self->layer_info[m].max_cumulative_bytes);
      extra_layers++;
      for (m=1; m < diff; m++)
        self->layer_info[n+m].max_cumulative_bytes = (int)
          exp((log_r1*(diff-m) + log_r2*m) / (double) diff);
    }
}

/*****************************************************************************/
/* STATIC                      create_quadtree                               */
/*****************************************************************************/

static void
  create_quadtree(int height, int width, code_subblock_ptr subblocks,
                  quad_siblings_ptr *first_entry,
                  quad_siblings_ptr *last_entry)

 /* Builds the quad-tree structure for a tree whose base measures
    `height' by `width'.  The base does not need to be square, nor do
    its dimensions need to be powers of 2.  Returns the head and tail of the
    array which represents the quad-tree via the `first_entry' and
    `last_entry' arguments. */

{
  int rows, cols, last_rows, last_cols, elts, leaves;
  int r, c, r_node, c_node, r_child, c_child, n, idx;
  quad_siblings_ptr head, tail, last_start, start, qp;

  rows = height; cols = width; elts = 0;
  while ((rows > 1) || (cols > 1))
    {
      rows = 1 + ((rows-1)>>1);
      cols = 1 + ((cols-1)>>1);
      elts += rows*cols;
    }
  if (elts == 0)
    {
      *first_entry = *last_entry = NULL;
      return;
    }
  head = (quad_siblings_ptr)
    local_malloc(sizeof(quad_siblings)*(size_t) elts);
  tail = head + (elts-1);
  last_start = tail + 1;
  last_rows = height; last_cols = width; leaves = 1;
  while ((last_rows > 1) || (last_cols > 1))
    {
      rows = 1 + ((last_rows-1)>>1);
      cols = 1 + ((last_cols-1)>>1);
      elts = rows*cols;
      start = last_start - elts;
      for (qp=start, r=0; r < rows; r++)
        for (c=0; c < cols; c++, qp++)
          {
            qp->flags = (leaves)?QUAD_LEAVES:0;
            qp->num_siblings = 0; /* Actual value will be set before coding. */
            for (n=0, r_node=0; r_node < 2; r_node++)
              for (c_node=0; c_node < 2; c_node++, n++)
                {
                  qp->nodes[n] = 0;
                  r_child = 2*r + r_node;
                  c_child = 2*c + c_node;
                  idx = r_child * last_cols + c_child;
                  if ((r_child >= last_rows) || (c_child >= last_cols))
                    qp->node_update[n] = NULL; /* Illegal node. */
                  else if (leaves)
                    qp->node_update[n] = (std_short *)(subblocks + idx);
                  else
                    qp->node_update[n] = (std_short *)(last_start + idx);
                }
          }
      leaves = 0;
      last_start = start;
      last_rows = rows;
      last_cols = cols;
    }
  assert(last_start == head);
  *first_entry = head;
  *last_entry = tail;
}

/* ========================================================================= */
/* ------------------------ Header Output Functions ------------------------ */
/* ========================================================================= */

/*****************************************************************************/
/* STATIC                   write_block_dimensions                           */
/*****************************************************************************/

static void
  write_block_dimensions(the_encoder_ref self)
{
  std_byte byte;
  int n;

  assert((self->master.max_dim <= 256) &&
         (self->master.max_dim >= self->master.max_sub_dim) &&
         (self->master.max_sub_dim > 0));
  byte = (std_byte)(self->master.max_dim & 0x00FF);
  self->output->push_bytes(self->output,&byte,1);
  byte = (std_byte)(self->master.max_sub_dim & 0x00FF);
  self->output->push_bytes(self->output,&byte,1);
  byte = (std_byte)(self->max_block_height & 0x00FF);
  self->output->push_bytes(self->output,&byte,1);
  for (n=0; n < self->bitstream_layers; n++)
    self->layer_info[n].max_cumulative_bytes -= 3;
}

/*****************************************************************************/
/* STATIC                      write_layer_info                              */
/*****************************************************************************/

static void
  write_layer_info(the_encoder_ref self)
{
  std_byte byte;
  int n;

  assert((self->bitstream_layers > 0) && (self->bitstream_layers <= 256));
  byte = (std_byte)(self->bitstream_layers);
  self->output->push_bytes(self->output,&byte,1);
  for (n=0; n < self->bitstream_layers; n++)
    self->layer_info[n].max_cumulative_bytes -= 1;
}


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

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

static void
  compute_rd_slopes(int num_points, double *cumulative_wmse_reduction,
                    pass_info_ptr passes, std_int *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);
          if (slope > ((double) (1<<30)))
            passes[n].rd_slope = 1<<30;
          else if (slope < 1.0)
            passes[n].rd_slope = 1;
          else
            passes[n].rd_slope = (std_int) floor(slope+0.5);
          if (n >= completed_points)
            completed_points = n+1;
        }
    }

  { /* Accumulate slope rates information. */
    std_int idx, threshold, slope_val, num_bytes;

    num_bytes = 0;
    for (n=0, idx=30, threshold=(1<<idx); idx >= 0; idx--, threshold>>=1)
      {
        for (; n < num_points; n++)
          {
            if ((slope_val = passes[n].rd_slope) >= threshold)
              num_bytes = passes[n].cumulative_bytes;
            else if (slope_val > 0)
              break;
          }
        slope_rates[idx] += num_bytes;
      }
  }
}

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

/*****************************************************************************/
/* STATIC                    transpose_sample_buffer                         */
/*****************************************************************************/

static void
  transpose_sample_buffer(ifc_int *src, ifc_int *dest,
                          int src_rows, int src_cols, int row_gap)

 /* Transposes the array in `src', copying the result to `dest'.  The source
    buffer has the indicated dimensions and for both buffers, the separation
    between rows is `row_gap', which must of course be no less than the
    maximum of `src_rows' and `src_cols'. */

{
  int r, c;
  ifc_int *sp, *dp, *spp, *dpp;

  assert((row_gap >= src_rows) && (row_gap >= src_cols));
  for (spp=src, dpp=dest, r=src_rows; r > 0; r--, spp+=row_gap, dpp++)
    for (sp=spp, dp=dpp, c=src_cols; c > 0; c--, sp++, dp+=row_gap)
      *dp = *sp;
}

/*****************************************************************************/
/* STATIC                      initialize_quadtree                           */
/*****************************************************************************/

static ifc_int
  initialize_quadtree(block_master_ptr master)

 /* Initializes the array of `quad_siblings' structures which runs from
    `first_quad' to `last_quad' in the `master' structure, to work with
    a `block_rows'x`block_cols' array of sample values.  The sample buffer
    itself is found in `master'->sample_buffer and has `master'->row_gap
    samples between consecutive rows of block samples.  Returns the
    inclusive OR of the magnitude bits from all samples in the block,
    even if the quad-tree is empty (i.e. `first_quad' and `last_quad' are
    both NULL). */

{
  quad_siblings_ptr scan;
  int num_siblings, b, rows, cols, r, c, row_gap;
  ifc_int val, *sp, *spp;

  row_gap = master->row_gap;
  if (master->first_quad == NULL)
    { /* Empty quad-tree case. */
      code_subblock_ptr subblock;

      subblock = master->first_subblock;
      assert(subblock == master->last_subblock);
      rows = subblock->rows;
      cols = subblock->cols;
      for (val=0, spp=master->sample_buffer, r=rows; r > 0; r--, spp+=row_gap)
        for (sp=spp, c=cols; c > 0; c--)
          val |= *(sp++);
      return(val & MAX_IFC_INT);
    }

  for (scan=master->last_quad; scan >= master->first_quad; scan--)
    if (scan->flags & QUAD_LEAVES)
      {
        scan->flags = QUAD_LEAVES; /* Reset all other flags. */
        for (num_siblings=b=0; b < 4; b++)
          {
            code_subblock_ptr subblock;

            subblock = (code_subblock_ptr)(scan->node_update[b]);
            if ((subblock == NULL) ||
                ((rows=subblock->rows) <= 0) ||
                ((cols=subblock->cols) <= 0))
              { scan->nodes[b] = MIN_IFC_INT; continue; }
            num_siblings++;
            for (val=0, spp=master->sample_buffer+subblock->offset,
                 r=rows; r > 0; r--, spp+=row_gap)
              for (sp=spp, c=cols; c > 0; c--)
                val |= *(sp++);
            scan->nodes[b] = val & MAX_IFC_INT;
          }
        scan->num_siblings = num_siblings;
      }
    else
      {
        scan->flags = 0; /* Reset flags. */
        for (num_siblings=b=0; b < 4; b++)
          {
            quad_siblings_ptr child;

            child = (quad_siblings_ptr)(scan->node_update[b]);
            if ((child == NULL) || (child->num_siblings == 0))
              { scan->nodes[b] = MIN_IFC_INT; continue; }
            num_siblings++;
            val = (child->nodes[0] | child->nodes[1] |
                   child->nodes[2] | child->nodes[3]) & MAX_IFC_INT;
            scan->nodes[b] = val;
          }
        scan->num_siblings = num_siblings;
      }
   scan = master->first_quad;
   val = (scan->nodes[0] | scan->nodes[1] |
          scan->nodes[2] | scan->nodes[3]) & MAX_IFC_INT;
   return(val);
}

/*****************************************************************************/
/* STATIC                        encode_quadtree                             */
/*****************************************************************************/

static void
  encode_quadtree(block_master_ptr master)

 /* Encodes the new quad-tree information for the current bit-plane (as
    identified by `master'->bit_idx). */

{
  quad_siblings_ptr scan, end;
  arith_state_ptr state;
  std_int A, C, word, rho;
  std_short avail, flags;
  ifc_int mask, val, symbol;
  int n;

  scan = master->first_quad;
  end = master->last_quad;
  if (scan == NULL)
    return; /* Nothing to do.  No quad-tree in use. */
  mask = 1 << master->bit_idx;
  state = &(master->coder_state);
  A=state->A; C=state->C; word=state->word; rho=state->rho;
  avail=state->available_bits;
  for (; scan <= end; scan++)
    {
      flags = scan->flags;
      if (!(flags & QUAD_PARENT_SIGNIFICANT))
        continue;
      if (!(flags & QUAD_NODE_SIGNIFICANT))
        { /* First time through. */
          int num_siblings;

          num_siblings = scan->num_siblings;
          for (n=0; n < 4; n++)
            if ((val=scan->nodes[n]) >= 0)
              { /* Otherwise, the node is not a valid sibling. */
                symbol = val & mask;
                num_siblings--;
                if ((num_siblings > 0) || (flags & QUAD_NODE_SIGNIFICANT))
                  { /* Significance of node must be coded. */
                    emit_uniform_symbol(state,A,C,word,avail,rho,symbol);
                  }
                else
                  { /* Node must be significant.  No need to code symbol. */
                    assert(symbol);
                  }
                if (symbol)
                  {
                    scan->nodes[n] = MIN_IFC_INT;
                    *(scan->node_update[n]) |= QUAD_PARENT_SIGNIFICANT;
                    flags |= QUAD_NODE_SIGNIFICANT;
                  }
              }
          assert(flags & QUAD_NODE_SIGNIFICANT);
          scan->flags = flags;
        }
      else
        { /* This is a refinement pass for this group of siblings. */
          for (n=0; n < 4; n++)
            if ((val=scan->nodes[n]) >= 0)
              { /* Otherwise, the node's significance is already known. */
                symbol = val & mask;
                emit_uniform_symbol(state,A,C,word,avail,rho,symbol);
                if (symbol)
                  {
                    scan->nodes[n] = MIN_IFC_INT;
                    *(scan->node_update[n]) |= QUAD_PARENT_SIGNIFICANT;
                  }
              }