lifting.c

JPEG2000实现的源码

/*****************************************************************************/
/* Author: Brendt Wohlberg (Los Alamos National Laboratory).                 */
/* Copyright 2001 University of California.                                  */
/*****************************************************************************/

#include <stdlib.h>
#include <assert.h>
#include <stdio.h>
#include "Lifting.h"
#include "Step.h"

extern int OverrideExtensionPolicy;
extern StdExtension EvenUpdateExtension;
extern StdExtension OddUpdateExtension;

int lift_fwd_irr(const FltBankIrr* fltbank, const float* input,
		 int input_offset, int input_length, float* output) {
  float *channels[2] = {0, 0};
  int channel_offsets[2] = {0, 0};
  int channel_lengths[2] = {0, 0};
  ExtensionDesc extensions[2] = {{{SExtType, WExtSymm, 1},
				  {SExtType, WExtSymm, 1}},
				 {{SExtType, WExtSymm, 1},
				  {SExtType, WExtSymm, 1}}};
  int s, m, k;

  /* Set extension policy */
  if (fltbank->Filt_Cat == 1) {/* WS filters */
    int left_boundary_odd, rght_boundary_odd, c;
    left_boundary_odd = input_offset%2;
    rght_boundary_odd = (input_offset+input_length-1)%2;
    for (c = 0; c < 2; c++) {
      extensions[c].leftext.extsymm = (left_boundary_odd^(1-c))?WExtSymm:
	HExtSymm;
      extensions[c].rghtext.extsymm = (rght_boundary_odd^(1-c))?WExtSymm:
	HExtSymm;
    }
  } else {
    if (OverrideExtensionPolicy) {
      assert(EvenUpdateExtension != NullStdExt &&
	     OddUpdateExtension != NullStdExt);
      extensions[0] = *stdext_to_extdesc(EvenUpdateExtension);
      extensions[1] = *stdext_to_extdesc(OddUpdateExtension);
    } else {
      extensions[0] = *stdext_to_extdesc(C00);
      extensions[1] = *stdext_to_extdesc(C00);
    }
  }

  channels[0] = malloc(sizeof(float) * input_length);
  if (!channels[0])
    return 0;

  /* Split input into polyphase channels */
  channel_offsets[0] = input_offset/2;
  channel_lengths[0] = input_length/2 + (input_length%2 && !(input_offset%2));
  channel_offsets[1] = (input_offset-1)/2;
  channel_lengths[1] = input_length/2 + (input_length%2 && input_offset%2);

  channels[1] = channels[0] + channel_lengths[0];
  for (k = 0; k < channel_lengths[0]; k++) {
    assert((2*k + input_offset%2) >= 0 && 
	   (2*k + input_offset%2) < input_length);
    channels[0][k] = input[2*k + input_offset%2];
  }
  for (k = 0; k < channel_lengths[1]; k++) {
    assert((2*k + 1-(input_offset%2)) >= 0 && 
	   (2*k + 1-(input_offset%2)) < input_length);
    channels[1][k] = input[2*k + 1-(input_offset%2)];
  }
  
  
  /* Work through lifting steps */
  for (s = 0, m = fltbank->m0; s < fltbank->N_LS; s++) {
    do_step_float(&extensions[1-m], fltbank->alpha[s], 
		  fltbank->off[s], fltbank->L[s],
		  channels[1-m], channel_offsets[1-m], channel_lengths[1-m],
		  channels[m], channel_offsets[m], channel_lengths[m]);
    m = 1-m;
  }

  /* Merge polyphase channels into output */
  for (k = 0; k < channel_lengths[0]; k++){
    assert((2*k + input_offset%2) >= 0 && 
	   (2*k + input_offset%2) < input_length);
    output[2*k + input_offset%2] = channels[0][k]/fltbank->K;
  }
  for (k = 0; k < channel_lengths[1]; k++) {
    assert((2*k + 1-(input_offset%2)) >= 0 && 
	   (2*k + 1-(input_offset%2)) < input_length);
    output[2*k + 1-(input_offset%2)] = fltbank->K*channels[1][k];
  }

  free(channels[0]);
  return 1;
}

int lift_inv_irr(const FltBankIrr* fltbank, const float* input,
		 int input_offset, int input_length, float* output) {
  float *channels[2] = {0, 0};
  int channel_offsets[2] = {0, 0};
  int channel_lengths[2] = {0, 0};
  ExtensionDesc extensions[2] = {{{SExtType, WExtSymm, 1},
				  {SExtType, WExtSymm, 1}},
				 {{SExtType, WExtSymm, 1},
				  {SExtType, WExtSymm, 1}}};
  int s, m, k;

  /* Set extension policy */
  if (fltbank->Filt_Cat == 1) {/* WS filters */
    int left_boundary_odd, rght_boundary_odd, c;
    left_boundary_odd = input_offset%2;
    rght_boundary_odd = (input_offset+input_length-1)%2;
    for (c = 0; c < 2; c++) {
      extensions[c].leftext.extsymm = (left_boundary_odd^(1-c))?WExtSymm:
	HExtSymm;
      extensions[c].rghtext.extsymm = (rght_boundary_odd^(1-c))?WExtSymm:
	HExtSymm;
    }
  } else {
    if (OverrideExtensionPolicy) {
      assert(EvenUpdateExtension != NullStdExt &&
	     OddUpdateExtension != NullStdExt);
      extensions[0] = *stdext_to_extdesc(EvenUpdateExtension);
      extensions[1] = *stdext_to_extdesc(OddUpdateExtension);
    } else {
      extensions[0] = *stdext_to_extdesc(C00);
      extensions[1] = *stdext_to_extdesc(C00);
    }
  }

  channels[0] = malloc(sizeof(float) * input_length);
  if (!channels[0])
    return 0;

  /* Split input into polyphase channels */
  channel_offsets[0] = input_offset/2;
  channel_lengths[0] = input_length/2 + (input_length%2 && !(input_offset%2));
  channel_offsets[1] = (input_offset-1)/2;
  channel_lengths[1] = input_length/2 + (input_length%2 && input_offset%2);  
  channels[1] = channels[0] + channel_lengths[0];
  for (k = 0; k < channel_lengths[0]; k++)
    channels[0][k] = fltbank->K*input[2*k + input_offset%2];
  for (k = 0; k < channel_lengths[1]; k++)
    channels[1][k] = input[2*k + 1-(input_offset%2)]/fltbank->K;
  
  
  /* Work through lifting steps */
  m = 1 - ((fltbank->m0 + fltbank->N_LS)%2);
  for (s = fltbank->N_LS-1; s >= 0; s--) {
    undo_step_float(&extensions[1-m], fltbank->alpha[s], 
		    fltbank->off[s], fltbank->L[s],
		    channels[1-m], channel_offsets[1-m], channel_lengths[1-m],
		    channels[m], channel_offsets[m], channel_lengths[m]);
    m = 1-m;
  }

  /* Merge polyphase channels into output */
  for (k = 0; k < channel_lengths[0]; k++)
    output[2*k + input_offset%2] = channels[0][k];
  for (k = 0; k < channel_lengths[1]; k++)
    output[2*k + 1-(input_offset%2)] = channels[1][k];

  free(channels[0]);
  return 1;
}


int lift_fwd_rev(const FltBankRev* fltbank, const int* input,
		 int input_offset, int input_length, int* output) {
  int *channels[2] = {0, 0};
  int channel_offsets[2] = {0, 0};
  int channel_lengths[2] = {0, 0};
  ExtensionDesc extensions[2] = {{{SExtType, WExtSymm, 1},
				  {SExtType, WExtSymm, 1}},
				 {{SExtType, WExtSymm, 1},
				  {SExtType, WExtSymm, 1}}};
  int s, m, k;

  /* Set extension policy */
  if (fltbank->Filt_Cat == 1) {/* WS filters */
    int left_boundary_odd, rght_boundary_odd, c;
    left_boundary_odd = input_offset%2;
    rght_boundary_odd = (input_offset+input_length-1)%2;
    for (c = 0; c < 2; c++) {
      extensions[c].leftext.extsymm = (left_boundary_odd^(1-c))?WExtSymm:
	HExtSymm;
      extensions[c].rghtext.extsymm = (rght_boundary_odd^(1-c))?WExtSymm:
	HExtSymm;
    }
  } else {
    if (OverrideExtensionPolicy) {
      assert(EvenUpdateExtension != NullStdExt &&
	     OddUpdateExtension != NullStdExt);
      extensions[0] = *stdext_to_extdesc(EvenUpdateExtension);
      extensions[1] = *stdext_to_extdesc(OddUpdateExtension);
    } else {
      extensions[0] = *stdext_to_extdesc(C00);
      extensions[1] = *stdext_to_extdesc(C00);
    }
  }

  channels[0] = malloc(sizeof(int) * input_length);
  if (!channels[0])
    return 0;

  /* Split input into polyphase channels */
  channel_offsets[0] = input_offset/2;
  channel_lengths[0] = input_length/2 + (input_length%2 && !(input_offset%2));
  channel_offsets[1] = (input_offset-1)/2;
  channel_lengths[1] = input_length/2 + (input_length%2 && input_offset%2);

  channels[1] = channels[0] + channel_lengths[0];
  for (k = 0; k < channel_lengths[0]; k++) {
    assert((2*k + input_offset%2) >= 0 && 
	   (2*k + input_offset%2) < input_length);
    channels[0][k] = input[2*k + input_offset%2];
  }
  for (k = 0; k < channel_lengths[1]; k++) {
    assert((2*k + 1-(input_offset%2)) >= 0 && 
	   (2*k + 1-(input_offset%2)) < input_length);
    channels[1][k] = input[2*k + 1-(input_offset%2)];
  }
  
  
  /* Work through lifting steps */
  for (s = 0, m = fltbank->m0; s < fltbank->N_LS; s++) {
    do_step_int(&extensions[1-m], fltbank->alpha[s],
		fltbank->off[s], fltbank->L[s],
		fltbank->epsilon[s], fltbank->beta[s],
		channels[1-m], channel_offsets[1-m], channel_lengths[1-m],
		channels[m], channel_offsets[m], channel_lengths[m]);
    m = 1-m;
  }

  /* Merge polyphase channels into output */
  for (k = 0; k < channel_lengths[0]; k++){
    assert((2*k + input_offset%2) >= 0 && 
	   (2*k + input_offset%2) < input_length);
    output[2*k + input_offset%2] = channels[0][k];
  }
  for (k = 0; k < channel_lengths[1]; k++) {
    assert((2*k + 1-(input_offset%2)) >= 0 && 
	   (2*k + 1-(input_offset%2)) < input_length);
    output[2*k + 1-(input_offset%2)] = channels[1][k];
  }

  free(channels[0]);
  return 1;
}

int lift_inv_rev(const FltBankRev* fltbank, const int* input,
		 int input_offset, int input_length, int* output) {
    int *channels[2] = {0, 0};
  int channel_offsets[2] = {0, 0};
  int channel_lengths[2] = {0, 0};
  ExtensionDesc extensions[2] = {{{SExtType, WExtSymm, 1},
				  {SExtType, WExtSymm, 1}},
				 {{SExtType, WExtSymm, 1},
				  {SExtType, WExtSymm, 1}}};
  int s, m, k;

  /* Set extension policy */
  if (fltbank->Filt_Cat == 1) {/* WS filters */
    int left_boundary_odd, rght_boundary_odd, c;
    left_boundary_odd = input_offset%2;
    rght_boundary_odd = (input_offset+input_length-1)%2;
    for (c = 0; c < 2; c++) {
      extensions[c].leftext.extsymm = (left_boundary_odd^(1-c))?WExtSymm:
	HExtSymm;
      extensions[c].rghtext.extsymm = (rght_boundary_odd^(1-c))?WExtSymm:
	HExtSymm;
    }
  } else {
    if (OverrideExtensionPolicy) {
      assert(EvenUpdateExtension != NullStdExt &&
	     OddUpdateExtension != NullStdExt);
      extensions[0] = *stdext_to_extdesc(EvenUpdateExtension);
      extensions[1] = *stdext_to_extdesc(OddUpdateExtension);
    } else {
      extensions[0] = *stdext_to_extdesc(C00);
      extensions[1] = *stdext_to_extdesc(C00);
    }
  }

  channels[0] = malloc(sizeof(int) * input_length);
  if (!channels[0])
    return 0;

  /* Split input into polyphase channels */
  channel_offsets[0] = input_offset/2;
  channel_lengths[0] = input_length/2 + (input_length%2 && !(input_offset%2));
  channel_offsets[1] = (input_offset-1)/2;
  channel_lengths[1] = input_length/2 + (input_length%2 && input_offset%2);  
  channels[1] = channels[0] + channel_lengths[0];
  for (k = 0; k < channel_lengths[0]; k++)
    channels[0][k] = input[2*k + input_offset%2];
  for (k = 0; k < channel_lengths[1]; k++)
    channels[1][k] = input[2*k + 1-(input_offset%2)];
  
  
  /* Work through lifting steps */
  m = 1 - ((fltbank->m0 + fltbank->N_LS)%2);
  for (s = fltbank->N_LS-1; s >= 0; s--) {
    undo_step_int(&extensions[1-m], fltbank->alpha[s], 
		  fltbank->off[s], fltbank->L[s],
		  fltbank->epsilon[s], fltbank->beta[s],
		  channels[1-m], channel_offsets[1-m], channel_lengths[1-m],
		  channels[m], channel_offsets[m], channel_lengths[m]);
    m = 1-m;
  }

  /* Merge polyphase channels into output */
  for (k = 0; k < channel_lengths[0]; k++)
    output[2*k + input_offset%2] = channels[0][k];
  for (k = 0; k < channel_lengths[1]; k++)
    output[2*k + 1-(input_offset%2)] = channels[1][k];

  free(channels[0]);
  return 1;
}