hplx_vert_analysis_by_lifting.c

JPEG2000实现的源码

      return(1);
    }
  ext_before = 0;
  if (min_row_idx < 0)
    {
      ext_before = -min_row_idx;
      num_rows -= ext_before;
      min_row_idx = 0;
    }
  ext_after = (num_rows+min_row_idx) - stage->rows;
  if (ext_after > 0)
    num_rows -= ext_after;
  else
    ext_after = 0;
  assert((num_rows > 0) && (min_row_idx >= stage->min_row_idx));
  step = stage->generator;
  while ((next_row_idx=stage->next_row_idx) < (min_row_idx+num_rows))
    { /* Need to generate this row. */
      hplx_row_buffer_state_ptr source_bufs[8], target_buf, output_buf;
      int source_rows, min_source_idx;
      int n, k;

      if (step == NULL)
        return(0); /* Need to wait until new rows are added from outside. */
      min_source_idx = next_row_idx - step->neg_support;
      min_source_idx += step->target->missing_first_row -
                        step->source->missing_first_row;
      source_rows = step->pos_support + step->neg_support + 1;
      assert(source_rows <= 8);
      if ((!hplx_lift_stage__get_rows(step->source,min_source_idx,source_rows,
                                      source_bufs)) ||
          (!hplx_lift_stage__get_rows(step->target,next_row_idx,1,
                                      &target_buf)))
        return(0);
      hplx_lift_stage__release_row(step->source,min_source_idx);
      hplx_lift_stage__release_row(step->target,next_row_idx);
      output_buf = hplx_lift_stage__add_row(stage);
      if (step->float_taps != NULL)
        { /* Apply floating-point lifting step. */
          float *taps, *sp[8], *tp, *op, sum;

          tp = target_buf->float_buf;
          op = output_buf->float_buf;
          if (source_bufs[0] == NULL)
            { /* Source range entirely outside image support; skip the step. */
              for (n=stage->cols; n > 0; n--)
                *(op++) = *(tp++);
            }
          else if (((next_row_idx == 0) && stage->special_first) ||
                   (next_row_idx == (stage->rows-stage->special_last)))
            { /* Simply scale the row instead of performing the actual
                 lifting step. */
              float dc_gain = step->dc_gain;

              for (n=stage->cols; n > 0; n--)
                *(op++) = *(tp++) * dc_gain;
            }
          else
            {
              for (k=0; k < source_rows; k++)
                sp[k] = source_bufs[k]->float_buf;
              taps = step->float_taps - step->neg_support;
              for (n=stage->cols; n > 0; n--)
                {
                  for (sum=0.0F, k=0; k < source_rows; k++)
                    sum += taps[k] * *((sp[k])++);
                  *(op++) = *(tp++) + sum;
                }
            }
        }
      else
        { /* Apply integer lifting step. */
          int *taps;
          ifc_int *sp[8], *tp, *op;
          std_int sum, offset, downshift;

          tp = target_buf->int_buf;
          op = output_buf->int_buf;
          if (source_bufs[0] == NULL)
            { /* Source range entirely outside image support; skip the step. */
              for (n=stage->cols; n > 0; n--)
                *(op++) = *(tp++);
            }
          else if (((next_row_idx == 0) && stage->special_first) ||
                   (next_row_idx == (stage->rows-stage->special_last)))
            { /* Skip lifting step on frame boundaries. */
              for (n=stage->cols; n > 0; n--)
                *(op++) = *(tp++);
            }
          else
            {
              for (k=0; k < source_rows; k++)
                sp[k] = source_bufs[k]->int_buf;
              taps = step->int_taps - step->neg_support;
              downshift = step->int_downshift;
              offset = (1<<downshift)>>1;
              for (n=stage->cols; n > 0; n--)
                {
                  for (sum=0, k=0; k < source_rows; k++)
                    sum += taps[k] * *((sp[k])++);
                  *(op++) = *(tp++) + ((sum+offset)>>downshift);
                }
            }
        }
    }

  /* Fill in the contents of the `bufs' array, applying the relevant
     symmetric extension policies as necessary. */

  offset = min_row_idx - stage->min_row_idx;
  fp = bufs + ext_before; /* Points to first interior row. */
  lp = fp + num_rows - 1; /* Points to last interior row. */
  for (n=0; n < num_rows; n++)
    fp[n] = stage->row_buffers[offset+n];

  /* Perform extensions, alternating between left and right extension
     until we are done.  Note that this code will probably not mimic
     symmetric extension for even length filters. */
  for (n=p=1; (n <= ext_before) || (p <= ext_after); )
    {
      int did_something;

      did_something = 0;
      if (n <= ext_before)
        { /* Still need extension to upper part of image. */
          if (stage->first_extension)
            { /* Half-point extension. */
              if (n <= (num_rows+p-1))
                { fp[-n] = fp[n-1]; n++; did_something = 1; }
            }
          else
            { /* Whole-point extension. */
              if (n < (num_rows+p-1))
                { fp[-n] = fp[n]; n++; did_something = 1; }
            }
        }
      if (p <= ext_after)
        { /* Still need extension to lower part of image. */
          if (stage->last_extension)
            { /* Half-point extension. */
              if (p <= (num_rows+n-1))
                { lp[p] = lp[-p+1]; p++; did_something = 1; }
            }
          else
            { /* Whole-point extension. */
              if (p < (num_rows+n-1))
                { lp[p] = lp[-p]; p++; did_something = 1; }
            }
        }
      assert(did_something);
    }
  return(1);
}


/* ========================================================================= */
/* ------------------ Implementation of Interface Functions ---------------- */
/* ========================================================================= */

/*****************************************************************************/
/* STATIC                         __initialize                               */
/*****************************************************************************/

static void
  __initialize(hplx_analysis_stage_ref base, forward_info_ref info,
               canvas_dims_ptr dims, canvas_dims_ptr frame_dims)
{
  hplx_vert_analysis_by_lifting_ref self =
    (hplx_vert_analysis_by_lifting_ref) base;
  frame_info frame;
  hplx_buffering_stage_ptr stage;
  hplx_lifting_step_ptr step;
  int neg_supports[INFO__MAX_LIFTING_STEPS];
  int pos_supports[INFO__MAX_LIFTING_STEPS];
  int int_downshifts[INFO__MAX_LIFTING_STEPS];
  int *int_taps[INFO__MAX_LIFTING_STEPS];
  float *float_taps[INFO__MAX_LIFTING_STEPS];
  int b, n;

  /* Begin by determining the Wavelet kernels and derivative properties. */

  self->kernel_type =
    info->get_kernel_type(info,base->direction,INFO__ANALYSIS,
                          base->component_idx,base->level_idx,NULL,NULL,
                          NULL,NULL);
  if (self->kernel_type == INFO__INT_LIFTING)
    {
      self->num_steps =
        info->get_int_steps(info,base->direction,base->component_idx,
                            base->level_idx,neg_supports,pos_supports,
                            int_taps,int_downshifts);
      for (b=0; b < base->vert_limit; b++)
        self->branch_needs_scaling[b] = 0;
    }
  else if (self->kernel_type == INFO__FLOAT_LIFTING)
    {
      self->num_steps =
        info->get_float_steps(info,base->direction,base->component_idx,
                              base->level_idx,neg_supports,pos_supports,
                              float_taps,self->branch_scale,
                              self->branch_scale+1);
      for (b=0; b < base->vert_limit; b++)
        {
          float val = self->branch_scale[b];
          self->branch_needs_scaling[b] = (val > 1.00001F) || (val < 0.99999F);
        }
      self->heap.use_floats = 1;
    }
  else
    assert(0);
  self->steps = (hplx_lifting_step_ptr)
    local_malloc(XFORM_MEM_KEY,sizeof(hplx_lifting_step)*self->num_steps);
  self->num_buffering_stages = self->num_steps + 2;
  self->buffering_stages = (hplx_buffering_stage_ptr)
    local_malloc(XFORM_MEM_KEY,
                 sizeof(hplx_buffering_stage)*self->num_buffering_stages);
  for (n=0; n < self->num_buffering_stages; n++)
    {
      stage = self->buffering_stages + n;
      stage->odd_lines = ((n & 1) == 0);
      stage->max_buffered_rows = 1;
      stage->num_users = 0;
      stage->row_buffers = (hplx_row_buffer_state_ptr *)
        local_malloc(XFORM_MEM_KEY,sizeof(hplx_row_buffer_state_ptr) *
                     stage->max_buffered_rows);
      if (n >= 2)
        stage->generator = self->steps + (n-2);
      stage->heap = &(self->heap);
    }
  for (n=0; n < self->num_steps; n++)
    {
      step = self->steps + n;
      step->source = self->buffering_stages + n + 1;
      step->target = self->buffering_stages + n;
      step->source->num_users++;
      step->target->num_users++;
      step->neg_support = neg_supports[n];
      step->pos_support = pos_supports[n];
      if (self->kernel_type == INFO__INT_LIFTING)
        {
          step->int_taps = int_taps[n];
          step->int_downshift = int_downshifts[n];
        }
      else
        step->float_taps = float_taps[n];
    }
  for (b=0; b < base->vert_limit; b++)
    {
      if (self->num_steps & 1)
        self->branch_output_stages[b] = self->buffering_stages +
          (self->num_buffering_stages - 2 + b);
      else
        self->branch_output_stages[b] = self->buffering_stages +
          (self->num_buffering_stages - 1 - b);
      assert(self->branch_output_stages[b]->odd_lines == b);
      self->branch_output_stages[b]->num_users++;
        /* The process which pulls lines out of the final stages is
           understood as a user of the relevant buffers. */
    }
  if (self->kernel_type == INFO__FLOAT_LIFTING)
    { /* Compute DC gains associated with the output of each update lifting
         step, i.e. the 1'st, 3'rd, ... lifting steps. */
      float current_gain, source_gain, step_gain;
      int k;

      for (n=0; n < self->num_steps; n++)
        {
          step = self->steps + n;
          current_gain = (n>=2)?(step[-2].dc_gain):1.0F;
          source_gain = (n>=1)?(step[-1].dc_gain):1.0F;
          step_gain = 0.0F;
          for (k=-step->neg_support; k <= step->pos_support; k++)
            step_gain += step->float_taps[k];
          step->dc_gain = current_gain + step_gain*source_gain;
        }
      for (current_gain=1.0F, n=1; n < self->num_steps; n+=2)
        { /* This step ensures that the update lifting steps identify
             the incremental DC gain relative to the output of the previous
             update step, rather than the absolute DC gain at the end of
             each step. */
          step = self->steps + n;
          step_gain = step->dc_gain / current_gain;
          current_gain = step->dc_gain;
          step->dc_gain = step_gain;
        }
    }

  /* Now determine tile and frame dimensions and attributes. */

  self->rows_left_in_tile = dims->rows;
  self->first_row_idx_in_frame = dims->top_row;
  self->heap.buffer_length = dims->cols;
  for (n=0; n < self->num_buffering_stages; n++)
    self->buffering_stages[n].cols = dims->cols;
  info->get_level_info(info,base->component_idx,base->level_idx,
                       NULL,NULL,&frame,NULL);
  if (frame.vert_ssodwt)
    {
      int uses_overlap;

      /* First check that filters use overlap. */
      for (uses_overlap=0, n=0; n < self->num_steps; n++)
        if (((n & 1)==0) && (self->steps[n].pos_support > 0))
          uses_overlap = 1; /* Prediction step uses overlapping boundary */
      if (!uses_overlap)
        local_error("The 1-sample overlap method has no effect whatsoever "
                    "if all update steps in the lifting implementation of "
                    "the filter have only two taps!  The method was designed "
                    "for use with odd length filters.");

      /* Now for the frame dimensions and location. */
      self->frame_rows = frame_dims->rows;
      self->frame_ref = frame_dims->top_row;
      self->first_frame_rows = frame_dims->rows -
        ((dims->top_row - frame_dims->top_row) & (frame_dims->rows - 1));
      if (self->first_frame_rows >= dims->rows)
        {
          self->first_frame_rows = dims->rows;
          self->num_frames = 1;
        }
      else
        self->num_frames = 2 +
          (dims->rows - self->first_frame_rows - 1) / frame_dims->rows;
    }
  else
    {
      self->num_frames = 1;
      self->frame_rows = self->first_frame_rows = dims->rows;
      self->frame_ref = dims->top_row;
    }

  /* Now perform branch initialization. */

  for (b=0; b < base->vert_limit; b++)
    {