ebcot_encoder.c

JPEG2000实现的源码

    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 -------------------------- */
/* ========================================================================= */

/*****************************************************************************/
/* STATIC                     apply_roi_boost_mask                           */
/*****************************************************************************/

static void
  apply_roi_boost_mask(ifc_int *sample_buffer, int rows, int cols,
                       int row_gap, int boost_delta, std_byte **boost_mask,
                       int mask_implicit)

 /* This function is called before coding a block for which the ROI
    boost values differ for different samples in the block.  The
    `boost_delta' argument identifies the difference between the minimum
    and maximum boost, while the `boost_mask' argument points to a
    two-dimensional array of values indicating the amount by which the actual
    boost for each sample exceeds the minimum boost.  The less significant
    samples have their magnitudes shifted down by `boost_delta' -
    `boost_mask[r][c]', where [r,c] is the index of the relevant sample
    within the code block.
        If `roi_implicit' is true (non-zero) we must zero out extra LSB's
    from the boosted regions to ensure that they will not be misdetected
    as background samples by the decoder.  This interferes slightly with
    the distortion estimation techniques, but it is the easiest way to
    handle the requirements of implicit ROI parititions. */

{
  ifc_int *sp, val, shift;
  std_byte *bp;
  int r, c;
  
  if (mask_implicit)
    {
      ifc_int lsb_mask;

      assert((boost_delta+1) <= IMPLEMENTATION_PRECISION);
      lsb_mask = (ifc_int)((-1) << (IMPLEMENTATION_PRECISION-(boost_delta+1)));
        /* Removes everything below the `boost_delta' most significant
           magnitude bits. */
      for (r=rows; r > 0; r--, boost_mask++, sample_buffer+=row_gap)
        for (sp=sample_buffer, bp=*boost_mask, c=cols; c > 0; c--, sp++, bp++)
          {
            shift = ((ifc_int) boost_delta) - ((ifc_int)(*bp));
            if (shift)
              {
                assert(shift == (ifc_int) boost_delta);
                val = *sp;
                val = (val & MIN_IFC_INT) | ((val & MAX_IFC_INT) >> shift);
                *sp = val;
              }
            else
              *sp &= lsb_mask;
          }
    }
  else
    {
      for (r=rows; r > 0; r--, boost_mask++, sample_buffer+=row_gap)
        for (sp=sample_buffer, bp=*boost_mask, c=cols; c > 0; c--, sp++, bp++)
          {
            shift = ((ifc_int) boost_delta) - ((ifc_int)(*bp));
            if (shift)
              {
                assert((shift >= 0) && (shift <= (ifc_int) boost_delta));
                val = *sp;
                val = (val & MIN_IFC_INT) | ((val & MAX_IFC_INT) >> shift);
                *sp = val;
              }
          }
    }
}

/*****************************************************************************/
/* STATIC                  interleave_sample_buffer                          */
/*****************************************************************************/

static ifc_int
  interleave_sample_buffer(block_master_ptr master, ifc_int *sample_buffer)

 /* Fills in the contents of the interleaved sample buffer from the
    supplied non-interleaved buffer, whose rows are separated by
    `master->sample_row_gap' array entries.  Returns
    the inclusive OR of all sample magnitudes. */

{
  int cols, rows, r, c, stripe_r, row_gap;
  ifc_int *sp, *dp, *spp, *dpp, val, block_mag;

  spp = sample_buffer;
  dpp = master->interleaved_sample_buffer;
  cols = master->width;
  rows = master->height;
  row_gap = master->sample_row_gap;
  block_mag = 0;
  for (stripe_r=4, r=rows; r > 0; r--, spp+=row_gap, dpp++, stripe_r--)
    {
      if (stripe_r == 0)
        { stripe_r=4; dpp += master->interleaved_row_gap - 4; }
      for (sp=spp, dp=dpp, c=cols; c > 0; c--, sp++, dp+=4)
        {
          val = *sp;
          block_mag |= val;
          *dp = val;
        }
    }
  block_mag &= MAX_IFC_INT;
  return(block_mag);
}

/*****************************************************************************/
/* STATIC                     reset_block_contexts                           */
/*****************************************************************************/

static void
  reset_block_contexts(block_master_ptr master)
{
  int rows, cols, c, d_gap;
  std_short *cp;
  std_int *dcp, *dcpp;
  std_int dbl_out_of_bounds; /* David T Cvis mod */

  rows = master->height;
  cols = master->width;
  dcpp = (std_int *)(master->interleaved_context_buffer);
  assert((master->interleaved_row_gap & 3) == 0);
  d_gap = master->interleaved_row_gap >> 1;
  dbl_out_of_bounds = (1<<OUT_OF_BOUNDS_POS);  /* David T Cvis mod */
  dbl_out_of_bounds <<= 16;                    /* David T Cvis mod */
  dbl_out_of_bounds |= (1<<OUT_OF_BOUNDS_POS); /* David T Cvis mod */
  for (; rows > 0; rows -= 4, dcpp+=d_gap)
    {
      for (dcp=dcpp, c=cols; c > 0; c--)
        { *(dcp++) = 0; *(dcp++) = 0; }
      /* Begin David T Cvis mod */
      for (c=8; c > 0; c--)
        { *(dcp++) = dbl_out_of_bounds; *(dcp++) = dbl_out_of_bounds; }
      /* End David T Cvis mod */
    }
  if (rows < 0)
    { /* One or more rows in the final stripe is out of bounds. */
      dcpp -= d_gap;
      for (; rows < 0; rows++)
        {
          cp = ((std_short *) dcpp) + (4+rows);
          for (c=cols; c > 0; c--, cp+=4)
            *cp = OUT_OF_BOUNDS;
        }
    }
}

/* Begin David T Cvis mod */
/*****************************************************************************/
/* STATIC                 compute_ll_band_local_activity                     */
/*****************************************************************************/

static void
  compute_ll_band_local_activity(ifc_int *src, ifc_int *dest, int rows,
                                 int cols, int row_gap)

 /* This function applies a high-pass filter to the samples in the `src'
    block to produce samples in the `dest' block, where the pass-band gain
    of the high-pass filter is equal to 1.0, so that the nominal dynamic
    range is unaffected by the operation.  Each sample in the `dest' array
    is set to half the magnitude of the difference between the corresponding
    sample in the `src' array and the average of its four immediate
    neighbours, bearing in mind that all values have a sign-magnitude
    representation.  The purpose of the function is to obtain samples
    for which the `compute_visibility_factors' function will produce
    meaningful results when the subband is the LL (i.e. DC) subband. */

{
  ifc_int *sp, *dp, top, bottom, left, right, centre, diff;
  int r, c;

  for (r=1; r <= rows; r++, src+=row_gap, dest+=row_gap)
    for (sp=src, dp=dest, c=1; c <= cols; c++, sp++, dp++)
      {
        centre = *sp;
        left = (c>1)?sp[-1]:centre;
        right =(c<cols)?sp[1]:centre;
        top = (r>1)?sp[-row_gap]:centre;
        bottom = (r<rows)?sp[row_gap]:centre;
        centre = (centre<0)?(MIN_IFC_INT-centre):centre;
        top = (top<0)?(MIN_IFC_INT-top):top;
        bottom = (bottom<0)?(MIN_IFC_INT-bottom):bottom;
        left = (left<0)?(MIN_IFC_INT-left):left;
        right = (right<0)?(MIN_IFC_INT-right):right;
        centre = (centre+1)>>1;
        top = (top+4)>>3;
        bottom = (bottom+4)>>3;
        left = (left+4)>>3;
        right = (right+4)>>3;
        diff = centre - (top+bottom+left+right);
        *dp = (diff<0)?-diff:diff;
      }
}
/* End David T Cvis mod */

/* Begin David T Cvis mod */
/*****************************************************************************/
/* STATIC                  compute_visibility_factors                        */
/*****************************************************************************/

static void
  compute_visibility_factors(block_master_ptr master, ifc_int *buffer,
                             double nominal_range,
                             double masked_sensitivity_reduction,
                             double max_masked_sensitivity_reduction)

 /* This function computes a visibility factor for each distortion cell which
    represents the amount by which distortion values for that cell should be
    scaled to account for their visibility in the context of masking by
    surrounding signals from the same subband.  The `nominal_range'
    argument identifies the nominal range of samples from the relevant
    subband, which is twice the expected maximum sample magnitude.
        The visibility model implemented by this function is as follows.
    The visual stimuli are assumed to be divided into visual bands which
    correspond to the individual subbands, b (not strictly true, but a useful
    approximation).  The stimulus in each subband is scaled by a CSF factor,
    C_b, which accounts mainly for the MTF of the optical and sensor
    integration processes in the HVS, as well as any assumed MTF for the
    display device (e.g. monitor or printer).  The sensitivity of the HVS
    to a sample, s_b[n], in subband b is then modeled as
                   (C_b)^2 * (w_b)^2 * (s_b[n])^2
        G * --------------------------------------------
            1 + (M_b-1) * (av_{k near n}(|s'_b[k]|^e))^2
    where G is a global factor, identical for all subbands, which is selected