www.pudn.com > EZW_.rar > FILTER.C
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
/* Port to C from C++
*
* Mow-Song, Ng 2/9/2002
* msng@mmu.edu.my
* http://www.pesona.mmu.edu.my/~msng
*
* I do not claim copyright to the code, but if you use them or modify them,
* please drop me a mail.
*
*/
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
/* Original copyright info */
/*---------------------------------------------------------------------------*/
// Baseline Wavelet Transform Coder Construction Kit
//
// Geoff Davis
// gdavis@cs.dartmouth.edu
// http://www.cs.dartmouth.edu/~gdavis
//
// Copyright 1996 Geoff Davis 9/11/96
//
// Permission is granted to use this software for research purposes as
// long as this notice stays attached to this software.
//
/*---------------------------------------------------------------------------*/
#include "filter.h"
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
/* Filter coefficients.
* The analysis and synthesis filters are grouped under FilterSet
*/
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
Real HaarCoeffs [] = { 1.0/SQRT2, 1.0/SQRT2 };
/* A few Daubechies filters */
Real Daub4Coeffs [] = { 0.4829629131445341, 0.8365163037378077,
0.2241438680420134, -0.1294095225512603 };
Real Daub6Coeffs [] = { 0.3326705529500825, 0.8068915093110924,
0.4598775021184914, -0.1350110200102546,
-0.0854412738820267, 0.0352262918857095 };
Real Daub8Coeffs [] = { 0.2303778133088964, 0.7148465705529154,
0.6308807679398587, -0.0279837694168599,
-0.1870348117190931, 0.0308413818355607,
0.0328830116668852, -0.0105974017850690 };
/* Filter from Eero Simoncelli's PhD thesis --
* used in Edward Adelson's EPIC wavelet coder
* These are probably the filter coefficients used in Shapiro's EZW paper
*/
Real AdelsonCoeffs[] = { 0.028220367, -0.060394127, -0.07388188,
0.41394752, 0.7984298, 0.41394752,
-0.07388188, -0.060394127, 0.028220367 };
/* 7/9 Filter from M. Antonini, M. Barlaud, P. Mathieu, and
* I. Daubechies, "Image coding using wavelet transform", IEEE
* Transactions on Image Processing", Vol. pp. 205-220, 1992.
*/
Real AntoniniSynthesis [] = { -6.453888262893856e-02,
-4.068941760955867e-02,
4.180922732222124e-01,
7.884856164056651e-01,
4.180922732222124e-01,
-4.068941760955867e-02,
-6.453888262893856e-02 };
Real AntoniniAnalysis[] = { 3.782845550699535e-02,
-2.384946501937986e-02,
-1.106244044184226e-01,
3.774028556126536e-01,
8.526986790094022e-01,
3.774028556126537e-01,
-1.106244044184226e-01,
-2.384946501937986e-02,
3.782845550699535e-02 };
// Unpublished 18/10 filter from Villasenor's group
Real Villa1810Synthesis [] = { 9.544158682436510e-04,
-2.727196296995984e-06,
-9.452462998353147e-03,
-2.528037293949898e-03,
3.083373438534281e-02,
-1.376513483818621e-02,
-8.566118833165798e-02,
1.633685405569902e-01,
6.233596410344172e-01,
6.233596410344158e-01,
1.633685405569888e-01,
-8.566118833165885e-02,
-1.376513483818652e-02,
3.083373438534267e-02,
-2.528037293949898e-03,
-9.452462998353147e-03,
-2.727196296995984e-06,
9.544158682436510e-04};
Real Villa1810Analysis [] = { 2.885256501123136e-02,
8.244478227504624e-05,
-1.575264469076351e-01,
7.679048884691438e-02,
7.589077294537618e-01,
7.589077294537619e-01,
7.679048884691436e-02,
-1.575264469076351e-01,
8.244478227504624e-05,
2.885256501123136e-02};
// Filters from Chris Brislawn's tutorial code
Real BrislawnAnalysis [] = { 0.037828455506995, -0.023849465019380,
-0.110624404418423, 0.377402855612654,
0.852698679009403,
0.377402855612654, -0.110624404418423,
-0.023849465019380, 0.037828455506995 };
Real BrislawnSynthesis [] = { -0.064538882628938, -0.040689417609558,
0.418092273222212,
0.788485616405664,
0.418092273222212,
-0.040689417609558, -0.064538882628938};
Real Brislawn2Analysis [] = { 0.026913419, -0.032303352,
-0.241109818, 0.054100420,
0.899506092, 0.899506092,
0.054100420, -0.241109818,
-0.032303352, 0.026913419};
Real Brislawn2Synthesis [] = { 0.019843545, 0.023817599,
-0.023257840, 0.145570740,
0.541132748, 0.541132748,
0.145570740, -0.023257840,
0.023817599, 0.019843545 };
/*
* Filters from J. Villasenor, B. Belzer, J. Liao, "Wavelet Filter
* Evaluation for Image Compression." IEEE Transactions on Image
* Processing, Vol. 2, pp. 1053-1060, August 1995.
*/
Real Villa1Analysis [] = { 3.782845550699535e-02,
-2.384946501937986e-02,
-1.106244044184226e-01,
3.774028556126536e-01,
8.526986790094022e-01,
3.774028556126537e-01,
-1.106244044184226e-01,
-2.384946501937986e-02,
3.782845550699535e-02};
Real Villa1Synthesis [] = { -6.453888262893856e-02,
-4.068941760955867e-02,
4.180922732222124e-01,
7.884856164056651e-01,
4.180922732222124e-01,
-4.068941760955867e-02,
-6.453888262893856e-02};
Real Villa2Analysis [] = { -8.472827741318157e-03,
3.759210316686883e-03,
4.728175282882753e-02,
-3.347508104780150e-02,
-6.887811419061032e-02,
3.832692613243884e-01,
7.672451593927493e-01,
3.832692613243889e-01,
-6.887811419061045e-02,
-3.347508104780156e-02,
4.728175282882753e-02,
3.759210316686883e-03,
-8.472827741318157e-03};
Real Villa2Synthesis [] = { 1.418215589126359e-02,
6.292315666859828e-03,
-1.087373652243805e-01,
-6.916271012030040e-02,
4.481085999263908e-01,
8.328475700934288e-01,
4.481085999263908e-01,
-6.916271012030040e-02,
-1.087373652243805e-01,
6.292315666859828e-03,
1.418215589126359e-02};
Real Villa3Analysis [] = { -1.290777652578771e-01,
4.769893003875977e-02,
7.884856164056651e-01,
7.884856164056651e-01,
4.769893003875977e-02,
-1.290777652578771e-01};
Real Villa3Synthesis [] = { 1.891422775349768e-02,
6.989495243807747e-03,
-6.723693471890128e-02,
1.333892255971154e-01,
6.150507673110278e-01,
6.150507673110278e-01,
1.333892255971154e-01,
-6.723693471890128e-02,
6.989495243807747e-03,
1.891422775349768e-02};
Real Villa4Analysis [] = { -1.767766952966369e-01,
3.535533905932738e-01,
1.060660171779821e+00,
3.535533905932738e-01,
-1.767766952966369e-01};
Real Villa4Synthesis [] = { 3.535533905932738e-01,
7.071067811865476e-01,
3.535533905932738e-01};
Real Villa5Analysis [] = { 7.071067811865476e-01,
7.071067811865476e-01};
Real Villa5Synthesis [] = { -8.838834764831845e-02,
8.838834764831845e-02,
7.071067811865476e-01,
7.071067811865476e-01,
8.838834764831845e-02,
-8.838834764831845e-02};
Real Villa6Analysis [] = { 3.314563036811943e-02,
-6.629126073623885e-02,
-1.767766952966369e-01,
4.198446513295127e-01,
9.943689110435828e-01,
4.198446513295127e-01,
-1.767766952966369e-01,
-6.629126073623885e-02,
3.314563036811943e-02};
Real Villa6Synthesis [] = { 3.535533905932738e-01,
7.071067811865476e-01,
3.535533905932738e-01};
Real OdegardAnalysis[] = { 5.2865768532960523e-02,
-3.3418473279346828e-02,
-9.3069263703582719e-02,
3.8697186387262039e-01,
7.8751377152779212e-01,
3.8697186387262039e-01,
-9.3069263703582719e-02,
-3.3418473279346828e-02,
5.2865768532960523e-02};
Real OdegardSynthesis[] = { -8.6748316131711606e-02,
-5.4836926902779436e-02,
4.4030170672498536e-01,
8.1678063499210640e-01,
4.4030170672498536e-01,
-5.4836926902779436e-02,
-8.6748316131711606e-02};
FILTERSET *Haar, *Daub4, *Daub6, *Daub8, *Antonini,
*Villa1810, *Adelson, *Brislawn, *Brislawn2,
*Villa1, *Villa2, *Villa3, *Villa4, *Villa5,
*Villa6, *Odegard;
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
/* FILTER */
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
FILTER * FilterAllocDef(void)
{
FILTER *filter;
/* allocate the filter */
if ((filter = (FILTER *)malloc(sizeof(FILTER)))==NULL){
return NULL;
}
filter->size=0;
filter->firstIndex=0;
filter->coeff=NULL;
return filter;
}
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
FILTER * FilterAlloc(int size, int firstIndex, Real *coeff)
{
FILTER *filter;
/* allocate the filter */
if ((filter = (FILTER *)malloc(sizeof(FILTER)))==NULL){
return NULL;
}
FilterInit(filter, size, firstIndex, coeff);
return filter;
}
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
FILTER * FilterAllocCopy(FILTER *filterSrc)
{
FILTER *filter;
/* allocate the filter */
if ((filter = (FILTER *)malloc(sizeof(FILTER)))==NULL){
return NULL;
}
filter->coeff=NULL;
FilterCopy(filter, filterSrc);
return filter;
}
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
void FilterDealloc(FILTER *filter)
{
if (filter->coeff!=NULL){
free(filter->coeff);
}
free(filter);
}
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
void FilterInit(FILTER *filter, int size, int firstIndex, Real *coeff)
{
int i;
filter->size=size;
filter->firstIndex=firstIndex;
filter->center=-firstIndex;
filter->coeff = (Real *)calloc(size, sizeof(Real));
if (coeff != NULL){
for (i=0; icoeff[i] = coeff[i];
}
}
else{
for (i=0; icoeff[i]=0;
}
}
}
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
Real FilterCoeff(FILTER *filter, int index)
{
return filter->coeff[index-filter->firstIndex];
}
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
void FilterCopy(FILTER *filterDest, FILTER *filterSrc)
{
if (filterDest->coeff!=NULL){
free(filterDest->coeff);
}
FilterInit(filterDest, filterSrc->size, filterSrc->firstIndex, filterSrc->coeff);
}
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
/* FILTERSET */
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
FILTERSET * FilterSetAllocDef(void)
{
FILTERSET *filterset;
if((filterset = (FILTERSET *)malloc(sizeof(FILTERSET)))==NULL){
return NULL;
}
filterset->symmetric=FALSE;
filterset->analysisLow=NULL;
filterset->analysisHigh=NULL;
filterset->synthesisLow=NULL;
filterset->synthesisHigh=NULL;
return filterset;
}
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
FILTERSET * FilterSetAlloc(int symmetric,
Real *anLow, int anLowSize, int anLowFirst,
Real *synLow, int synLowSize, int synLowFirst)
{
FILTERSET *filterset;
int i, sign;
if((filterset = (FILTERSET *)malloc(sizeof(FILTERSET)))==NULL){
return NULL;
}
/* 7/9/2002 */
filterset->symmetric = symmetric;
filterset->analysisLow = FilterAlloc(anLowSize, anLowFirst, anLow);
//- TO DO : Error checking
// If no synthesis coeffs are given, assume wavelet is orthogonal
if (synLow == NULL){
filterset->synthesisLow=FilterAllocCopy(filterset->analysisLow);
//- TO DO : Error checking
// For orthogonal wavelets, compute the high pass filter using
// the relation g_n = (-1)^n h_{1-n}^*
// (or equivalently g_{1-n} = (-1)^{1-n} h_n^*)
filterset->analysisHigh = FilterAlloc(filterset->analysisLow->size,
2 - filterset->analysisLow->size -
filterset->analysisLow->firstIndex, NULL);
//- TO DO : Error checking
// Compute (-1)^(1-n) for first n
if (filterset->analysisLow->firstIndex % 2){
sign = 1;
}
else{
sign = -1;
}
for (i = 0; i < filterset->analysisLow->size; i++){
filterset->analysisHigh->coeff[1 - i - filterset->analysisLow->firstIndex -
filterset->analysisHigh->firstIndex] = sign * filterset->analysisLow->coeff[i];
assert(1-i-filterset->analysisLow->firstIndex-filterset->analysisHigh->firstIndex >= 0);
assert (1-i-filterset->analysisLow->firstIndex - filterset->analysisHigh->firstIndex
< filterset->analysisHigh->size);
sign *= -1;
}
// Copy the high pass analysis filter to the synthesis filter
filterset->synthesisHigh = FilterAllocCopy(filterset->analysisHigh);
//- TO DO : Error checking
}
else{
// If separate synthesis coeffs given, assume biorthogonal
filterset->synthesisLow = FilterAlloc(synLowSize, synLowFirst, synLow);
//- TO DO : Error checking
// For orthogonal wavelets, compute the high frequency filter using
// the relation g_n = (-1)^n complement (h~_{1-n}) and
// g~_n = (-1)^n complement (h_{1-n})
// (or equivalently g_{1-n} = (-1)^{1-n} complement (h~_n))
filterset->analysisHigh = FilterAlloc(filterset->synthesisLow->size,
2 - filterset->synthesisLow->size - filterset->synthesisLow->firstIndex, NULL);
//- TO DO : Error checking
// Compute (-1)^(1-n) for first n
if (filterset->synthesisLow->firstIndex % 2){
sign = 1;
}
else{
sign = -1;
}
for (i = 0; i < filterset->synthesisLow->size; i++){
filterset->analysisHigh->coeff[1 - i - filterset->synthesisLow->firstIndex -
filterset->analysisHigh->firstIndex] = sign * filterset->synthesisLow->coeff[i];
assert (1 - i - filterset->synthesisLow->firstIndex -
filterset->analysisHigh->firstIndex >= 0);
assert (1 - i - filterset->synthesisLow->firstIndex -
filterset->analysisHigh->firstIndex < filterset->analysisHigh->size);
sign *= -1;
}
filterset->synthesisHigh = FilterAlloc(filterset->analysisLow->size,
2 - filterset->analysisLow->size - filterset->analysisLow->firstIndex, NULL);
//- TO DO : Error checking
// Compute (-1)^(1-n) for first n
if (filterset->analysisLow->firstIndex % 2){
sign = 1;
}
else{
sign = -1;
}
for (i = 0; i < filterset->analysisLow->size; i++) {
filterset->synthesisHigh->coeff[1 - i - filterset->analysisLow->firstIndex -
filterset->synthesisHigh->firstIndex] = sign * filterset->analysisLow->coeff[i];
assert (1 - i - filterset->analysisLow->firstIndex -
filterset->synthesisHigh->firstIndex >= 0);
assert (1 - i - filterset->analysisLow->firstIndex -
filterset->synthesisHigh->firstIndex < filterset->synthesisHigh->size);
sign *= -1;
}
}
return filterset;
}
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
FILTERSET * FilterSetAllocCopy(FILTERSET *filtersetSrc)
{
FILTERSET *filterset;
if((filterset = (FILTERSET *)malloc(sizeof(FILTERSET)))==NULL){
return NULL;
}
FilterSetCopy(filterset, filtersetSrc);
return filterset;
}
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
void FilterSetDealloc(FILTERSET *filterset)
{
FilterDealloc(filterset->analysisLow);
FilterDealloc(filterset->analysisHigh);
FilterDealloc(filterset->synthesisLow);
FilterDealloc(filterset->synthesisHigh);
free(filterset);
}
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
void FilterSetCopy(FILTERSET *filtersetDest, FILTERSET *filtersetSrc)
{
filtersetDest->symmetric=filtersetSrc->symmetric;
filtersetDest->analysisLow = FilterAllocCopy(filtersetSrc->analysisLow);
filtersetDest->analysisHigh = FilterAllocCopy(filtersetSrc->analysisHigh);
filtersetDest->synthesisLow = FilterAllocCopy(filtersetSrc->synthesisLow);
filtersetDest->synthesisHigh = FilterAllocCopy(filtersetSrc->synthesisHigh);
}
/*----------------------------------------------------------------------------*/
/*----------------------------------------------------------------------------*/
void InitializeFilterSets(void)
{
Haar = FilterSetAlloc(FALSE, HaarCoeffs, 2, 0, NULL, 0, 0);
Daub4 = FilterSetAlloc(FALSE, Daub4Coeffs, 4, 0, NULL, 0, 0);
Daub6 = FilterSetAlloc(FALSE, Daub6Coeffs, 6, 0, NULL, 0, 0);
Daub8 = FilterSetAlloc(FALSE, Daub8Coeffs, 8, 0, NULL, 0, 0);
Antonini = FilterSetAlloc(TRUE, AntoniniAnalysis, 9, -4,
AntoniniSynthesis, 7, -3);
Villa1810 = FilterSetAlloc(TRUE, Villa1810Analysis, 10, -4,
Villa1810Synthesis, 18, -8);
Adelson = FilterSetAlloc(TRUE, AdelsonCoeffs, 9, -4, NULL, 0, 0);
Brislawn = FilterSetAlloc(TRUE, BrislawnAnalysis, 9, -4,
BrislawnSynthesis, 7, -3);
Brislawn2 = FilterSetAlloc(TRUE, Brislawn2Analysis, 10, -4,
Brislawn2Synthesis, 10, -4);
Villa1 = FilterSetAlloc(TRUE, Villa1Analysis, 9, -4, Villa1Synthesis, 7, -3);
Villa2 = FilterSetAlloc(TRUE, Villa2Analysis, 13, -6, Villa2Synthesis, 11, -5);
Villa3 = FilterSetAlloc(TRUE, Villa3Analysis, 6, -2, Villa3Synthesis, 10, -4);
Villa4 = FilterSetAlloc(TRUE, Villa4Analysis, 5, -2, Villa4Synthesis, 3, -1);
Villa5 = FilterSetAlloc(TRUE, Villa5Analysis, 2, 0, Villa5Synthesis, 6, -2);
Villa6 = FilterSetAlloc(TRUE, Villa6Analysis, 9, -4, Villa6Synthesis, 3, -1);
Odegard = FilterSetAlloc(TRUE, OdegardAnalysis, 9, -4, OdegardSynthesis, 7, -3);
}
/*----------------------------------------------------------------------------*/
/*----------------------------------------------------------------------------*/
void RemoveFilterSets(void)
{
FilterSetDealloc(Haar);
FilterSetDealloc(Daub4);
FilterSetDealloc(Daub6);
FilterSetDealloc(Daub8);
FilterSetDealloc(Antonini);
FilterSetDealloc(Villa1810);
FilterSetDealloc(Adelson);
FilterSetDealloc(Brislawn);
FilterSetDealloc(Brislawn2);
FilterSetDealloc(Villa1);
FilterSetDealloc(Villa2);
FilterSetDealloc(Villa3);
FilterSetDealloc(Villa4);
FilterSetDealloc(Villa5);
FilterSetDealloc(Villa6);
FilterSetDealloc(Odegard);
}
/*----------------------------------------------------------------------------*/
/*----------------------------------------------------------------------------*/
void FilterWarning(char *fmt, ...)
{
va_list argptr;
va_start( argptr, fmt );
fprintf( stderr, "FilterWarning: " );
vprintf( fmt, argptr );
va_end( argptr );
}
/*----------------------------------------------------------------------------*/
/*----------------------------------------------------------------------------*/
void FilterError(char *fmt, ...)
{
va_list argptr;
va_start( argptr, fmt );
fprintf(stderr, "FilterError: " );
vprintf( fmt, argptr );
va_end( argptr );
exit( -1 );
}