sfq.c

spiht for linux this is used to decod and encode vedio i wich all enjoy

/*
 * 
 * QccPack: Quantization, compression, and coding libraries
 * Copyright (C) 1997-2009  James E. Fowler
 * 
 * This library is free software; you can redistribute it and/or
 * modify it under the terms of the GNU Library General Public
 * License as published by the Free Software Foundation; either
 * version 2 of the License, or (at your option) any later version.
 * 
 * This library is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * Library General Public License for more details.
 * 
 * You should have received a copy of the GNU Library General Public
 * License along with this library; if not, write to the
 * Free Software Foundation, Inc., 675 Mass Ave, Cambridge,
 * MA 02139, USA.
 * 
 */


#include "libQccPack.h"


int QccWAVsfqWaveletAnalysis(const QccIMGImageComponent *input_image,
                             QccWAVSubbandPyramid *subband_pyramid,
                             QccWAVZerotree *zerotree,
                             const QccSQScalarQuantizer *mean_quantizer,
                             const QccWAVWavelet *wavelet,
                             const QccWAVPerceptualWeights *perceptual_weights,
                             double lambda)
{
  int return_value;
  int num_levels;

  if ((subband_pyramid == NULL) ||
      (mean_quantizer == NULL) ||
      (wavelet == NULL))
    return(0);

  num_levels = zerotree->num_levels;

  if (QccMatrixCopy(subband_pyramid->matrix,
                    input_image->image,
                    input_image->num_rows,
                    input_image->num_cols))
    {
      QccErrorAddMessage("(QccWAVsfqWaveletAnalysis): Error calling QccMatrixCopy()");
      goto Error;
    }
  
  if (QccWAVSubbandPyramidSubtractMean(subband_pyramid,
                                       &(zerotree->image_mean),
                                       mean_quantizer))
    {
      QccErrorAddMessage("(QccWAVsfqWaveletAnalysis): Error calling QccWAVSubbandPyramidSubtractMean()");
      goto Error;
    }
  
  if (num_levels)
    if (QccWAVSubbandPyramidDWT(subband_pyramid,
                                num_levels,
                                wavelet))
      {
        QccErrorAddMessage("(QccWAVsfqWaveletAnalysis): Error calling QccWAVSubbandPyramidDWT()");
        goto Error;
      }

  if (perceptual_weights != NULL)
    if (QccWAVPerceptualWeightsApply(subband_pyramid,
                                     perceptual_weights))
      {
        QccErrorAddMessage("(QccWAVsfqWaveletAnalysis): Error calling QccWAVPerceptualWeightsApply()");
        QccErrorExit();
      }

  return_value = 0;
  goto Return;
 Error:
  return_value = 1;
 Return:
  return(return_value);
}


static int QccWAVsfqBasebandEncodeProcess(QccIMGImageComponent *baseband_image,
                                          QccSQScalarQuantizer 
                                          *baseband_quantizer,
                                          QccChannel *channels,
                                          QccIMGImageComponent 
                                          *reconstructed_baseband,
                                          double *distortion, 
                                          double *rate)
{
  int return_value;
  QccVector distortion_vector = NULL;

  if ((baseband_image == NULL) ||
      (baseband_quantizer == NULL) ||
      (channels == NULL))
    return(0);

  if ((distortion_vector = QccVectorAlloc(baseband_image->num_rows *
                                          baseband_image->num_cols)) == 
      NULL)
    {
      QccErrorAddMessage("(QccWAVsfqBasebandEncodeProcess): Error calling QccVectorAlloc()");
      goto Error;
    }

  if (QccIMGImageComponentScalarQuantize(baseband_image,
                                         baseband_quantizer,
                                         distortion_vector,
                                         &(channels[0])))
    {
      QccErrorAddMessage("(QccWAVsfqBasebandEncodeProcess): Error calling QccIMGImageComponentScalarQuantize()");
      goto Error;
    }

  if (QccIMGImageComponentInverseScalarQuantize(&(channels[0]),
                                                baseband_quantizer,
                                                reconstructed_baseband))
    {
      QccErrorAddMessage("(QccWAVsfqBasebandEncodeProcess): Error calling QccIMGImageComponentInverseScalarQuantize()");
      goto Error;
    }

  if (QccChannelNormalize(&(channels[0])))
    {
      QccErrorAddMessage("(QccWAVsfqBasebandEncodeProcess): Error calling QccChannelNormalize()");
      goto Error;
    }

  if (distortion != NULL)
    *distortion = 
      QccVectorMean(distortion_vector,
                    baseband_image->num_rows *
                    baseband_image->num_cols);
  if (rate != NULL)
    *rate = QccChannelEntropy(&(channels[0]), 1);

  return_value = 0;
  goto Return;
 Error:
  return_value = 1;
 Return:
  QccVectorFree(distortion_vector);
  return(return_value);
}


int QccWAVsfqBasebandEncode(const QccWAVSubbandPyramid *subband_pyramid,
                            QccWAVZerotree *zerotree,
                            QccSQScalarQuantizer *baseband_quantizer,
                            QccChannel *channels,
                            double lambda,
                            QccIMGImageComponent *reconstructed_baseband)
{
  int return_value;
  QccIMGImageComponent baseband_image;
  double distortion = 0;
  double rate = 0;
  double baseband_quantizer_start_stepsize;
  double baseband_quantizer_end_stepsize;
  int done = 0;
  int last_time = 0;
  QccVector J = NULL;
  int pass;
  int num_passes;
  double min_J = MAXDOUBLE;
  double winner = 0;
  double max_coefficient = -MAXDOUBLE;
  int row, col;

  if ((subband_pyramid == NULL) ||
      (zerotree == NULL) ||
      (baseband_quantizer == NULL) ||
      (channels == NULL))
    return(0);

  if (baseband_quantizer->stepsize <= 0.0)
    {
      baseband_quantizer_start_stepsize =
        QCCWAVSFQ_BASEBANDQUANTIZER_START_STEPSIZE;
      baseband_quantizer_end_stepsize =
        QCCWAVSFQ_BASEBANDQUANTIZER_END_STEPSIZE;
    }
  else
    {
      baseband_quantizer_start_stepsize =
        baseband_quantizer_end_stepsize =
        baseband_quantizer->stepsize;
      last_time = 1;
    }
  
  num_passes =
    ceil((baseband_quantizer_end_stepsize -
          baseband_quantizer_start_stepsize) /
         QCCWAVSFQ_BASEBANDQUANTIZER_INCREMENT) + 1;

  baseband_image.num_rows = zerotree->num_rows[0];
  baseband_image.num_cols = zerotree->num_cols[0];
  baseband_image.image = subband_pyramid->matrix;
  for (row = 0; row < baseband_image.num_rows; row++)
    for (col = 0; col < baseband_image.num_cols; col++)
      if (fabs(baseband_image.image[row][col]) > max_coefficient)
        max_coefficient = fabs(baseband_image.image[row][col]);

  if ((J = QccVectorAlloc(num_passes)) == NULL)
    {
      QccErrorAddMessage("(QccWAVsfqBasebandEncode): Error calling QccVectorAlloc()");
      goto Error;
    }
  for (pass = 0; pass < num_passes; pass++)
    J[pass] = MAXDOUBLE;

  baseband_quantizer->stepsize = 
    baseband_quantizer_start_stepsize;
  
  pass = 0;
  do
    {
      if (last_time)
        done = 1;
      
      if (QccSQScalarQuantization(max_coefficient,
                                  baseband_quantizer,
                                  NULL,
                                  &(channels[0].alphabet_size)))
        {
          QccErrorAddMessage("(QccWAVsfqBasebandEncode): Error calling QccSQScalarQuantization()");
          goto Error;
        }
      channels[0].alphabet_size =
        channels[0].alphabet_size * 2 + 1;
      baseband_quantizer->num_levels = channels[0].alphabet_size;

      if (QccWAVsfqBasebandEncodeProcess(&baseband_image,
                                         baseband_quantizer,
                                         channels,
                                         reconstructed_baseband,
                                         &distortion, &rate))
        {
          QccErrorAddMessage("(QccWAVsfqBasebandEncode): Error calling QccWAVsfqBasebandEncodeProcess()");
          goto Error;
        }

      if (!done)
        {
          J[pass] =
            distortion + lambda * rate;

          /****/
          /*
            printf("    step: %f\n", baseband_quantizer->stepsize);
            printf("     MSE: %f\n", distortion);
            printf("    rate: %f\n", rate);
            printf("       J: %f\n\n", J[pass]);
          */

          if (pass == num_passes - 1)
            {
              for (pass = 0;
                   pass < num_passes;
                   pass++)
                if (J[pass] < min_J)
                  {
                    min_J = J[pass];
                    winner = baseband_quantizer_start_stepsize +
                      QCCWAVSFQ_BASEBANDQUANTIZER_INCREMENT * pass;
                  }
              
              baseband_quantizer->stepsize = winner;
              last_time = 1;
            }
          else
            {
              baseband_quantizer->stepsize +=
                QCCWAVSFQ_BASEBANDQUANTIZER_INCREMENT;
              pass++;
            }
        }
    }
  while (!done);

  return_value = 0;
  goto Return;
 Error:
  return_value = 1;
 Return:
  QccVectorFree(J);
  return(return_value);
}


static int QccWAVsfqCalcCodewordLengths(QccVector codeword_lengths,
                                        int num_symbols,
                                        QccChannel *channels,
                                        QccWAVZerotree *zerotree)
{
  int subband, row, col;
  int symbol;
  int channel_index;
  int num_nonnull_symbols = 0;
  
  for (symbol = 0; symbol < num_symbols; symbol++)
    codeword_lengths[symbol] = 0;
  
  for (subband = 1;
       subband < zerotree->num_subbands; subband++)
    for (row = 0, channel_index = 0; 
         row < zerotree->num_rows[subband]; row++)
      for (col = 0; col < zerotree->num_cols[subband]; 
           col++, channel_index++)
        if (!QccWAVZerotreeNullSymbol(zerotree->zerotree[subband][row][col]))
          {
            codeword_lengths[channels[subband].channel_symbols
                            [channel_index]]++;
            num_nonnull_symbols++;
          }
  
  for (symbol = 0; symbol < num_symbols; symbol++)
    codeword_lengths[symbol] =
      -QccMathLog2(codeword_lengths[symbol] / num_nonnull_symbols);
  
  return(0);
}


static double QccWAVsfqCalcChildrenCost(QccVector *costs,
                                        int child_subband,
                                        int parent_row,
                                        int parent_col,
                                        int num_child_cols)
{
  int child_row, child_col;
  int num_children;
  double cost = 0;
  
  num_children = (child_subband < 4) ? 1 : 2;
  
  for (child_row = num_children*parent_row;
       child_row < num_children*(parent_row + 1);
       child_row++)
    for (child_col = num_children*parent_col;
         child_col < num_children*(parent_col + 1);
         child_col++)
      cost +=
        costs[child_subband][child_row*num_child_cols + child_col];
  
  return(cost);
  
}

static double QccWAVsfqSelectMinimumCost(QccWAVZerotree *zerotree, 
                                         int subband, 
                                         int row, int col, int cost_index,
                                         QccVector *J, 
                                         QccVector *residue_tree_J,
                                         QccVector *delta_J,
                                         QccVector *squared_subband_images,
                                         int *zerotree_pruned)
{
  int level;
  double cost1, cost2;
  double minimum_cost;
  int child_subband;
  int child_subband_start, child_subband_end;
  
  level = QccWAVSubbandPyramidCalcLevelFromSubband(subband, zerotree->num_levels);
  
  if (subband)
    child_subband_start = child_subband_end =
      subband + 3;
  else