#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <string.h>
#include <time.h>
#define F_EXTPAD 4
#define D_EXTPAD 2
/* float 9/7
you can get the coefficient by using [a,b,c,d]=WFILTERS('bior4.4') in matlab
********************************/
double h[5] ={ 0.85269867900889 , 0.37740285561283 , -0.11062440441844, -0.02384946501956 , 0.03782845550726 };
double g[4] ={ -0.78848561640558 , 0.41809227322162 , 0.04068941760916, -0.06453888262870 };
double h_b[4] ={ 0.78848561640558 , 0.41809227322162 , -0.04068941760916, -0.06453888262870 };
double g_b[5] ={-0.85269867900889 , 0.37740285561283 , 0.11062440441844, -0.02384946501956 , -0.03782845550726 };
/********* float 5/3 **********
you can get the coefficient by using [a,b,c,d]=WFILTERS('bior2.2') in matlab
**********************************
double h[5] ={ 1.06066017177982 , 0.35355339059327 ,-0.17677669529664 ,0};
double g[4] ={ -0.70710678118655 , 0.35355339059327 , 0 };
double h_b[4] ={ 0.70710678118655 , 0.35355339059327, 0 };
double g_b[5] ={-1.06066017177982 , 0.35355339059327, 0.17677669529664,0 };
*/
static float *x_alloc = NULL; // **
unsigned char **Alloc_Char_Matrix(short r, short c)
{
unsigned char *x, **y;
short n;
// calloc
x = (unsigned char *) calloc((long) r * c, sizeof(unsigned char));
y = (unsigned char **) calloc(r, sizeof(unsigned char *));
for (n = 0; n <= r - 1; n++)
y[n] = &x[(long) c * n];
return y;
}
short **Alloc_Short_Matrix(short r, short c)
{
short *x, **y;
short n;
// calloc
x = (short *) calloc((long) r * c, sizeof(short));
y = (short **) calloc(r, sizeof(short *));
for (n = 0; n <= r - 1; n++)
y[n] = &x[(long) c * n];
return y;
}
float **Alloc_Float_Matrix(short r, short c)
{
float *x, **y;
short n;
// calloc
x = (float *) calloc((long) r * c, sizeof(float));
y = (float **) calloc(r, sizeof(float *));
for (n = 0; n <= r - 1; n++)
y[n] = &x[(long) c * n];
return y;
}
void File2_Char_Matrix(unsigned char *f_name, unsigned char **x, short r, short c)
{
short i, j;
FILE *f_in;
f_in = fopen(f_name, "rb");
for (i = 0; i < r; i++)
for (j = 0; j < c; j++)
fread(&x[i][j], sizeof(unsigned char), 1, f_in);
fclose(f_in);
}
void Char_Matrix_2File(unsigned char *f_name, unsigned char **x, short r,short c)
{
short i, j;
FILE *f_out;
f_out = fopen(f_name, "wb");
for (i = 0; i < r; i++)
for (j = 0; j < c; j++)
fwrite(&x[i][j], sizeof(unsigned char), 1, f_out);
fclose(f_out);
}
void File2_Short_Matrix(unsigned char *f_name, short **x, short r, short c)
{
short i, j;
FILE *f_in;
f_in = fopen(f_name, "rb");
for (i = 0; i < r; i++)
for (j = 0; j < c; j++)
fread(&x[i][j], sizeof(short), 1, f_in);
fclose(f_in);
}
void Short_Matrix_2File(unsigned char *f_name, short **x, short r,short c)
{
short i, j;
FILE *f_out;
f_out = fopen(f_name, "wb");
for (i = 0; i < r; i++)
for (j = 0; j < c; j++)
fwrite(&x[i][j], sizeof(short), 1, f_out);
fclose(f_out);
}
void File2_Float_Matrix(unsigned char *f_name, float **x, short r, short c)
{
short i, j;
FILE *f_in;
f_in = fopen(f_name, "rb");
for (i = 0; i < r; i++)
for (j = 0; j < c; j++)
fread(&x[i][j], sizeof(float), 1, f_in);
fclose(f_in);
}
void Float_Matrix_2File(unsigned char *f_name, float **x, short r,short c)
{
short i, j;
FILE *f_out;
f_out = fopen(f_name, "wb");
for (i = 0; i < r; i++)
for (j = 0; j < c; j++)
fwrite(&x[i][j], sizeof(float), 1, f_out);
fclose(f_out);
}
void forwardf97(float *x_in, short N)
{
short i, n, half;
float *x, *r, *d;
x = x_alloc + F_EXTPAD;
memcpy(x, x_in, sizeof(float) * N);
for (i = 1; i <= F_EXTPAD; i++)
{
x[-i] = x[i];
x[(N - 1) + i] = x[(N - 1) - i];
}
half = N >> 1;
r = x_in;
d = x_in + half;
for (n = half; n--;)
{
*r++ = h[4] * (x[4] + x[-4]) + h[3] * (x[3] + x[-3]) + h[2] * (x[2] + x[-2]) +
h[1] * (x[1] + x[-1]) + h[0] * x[0];
x++;
*d++ = g[3] * (x[3] + x[-3]) + g[2] * (x[2] + x[-2]) + g[1] * (x[1] + x[-1]) +
g[0] * x[0];
x++;
}
}
void inversef97(float *x, short N)
{
short i, n, half;
float *r, *d;
half = N / 2;
r = x_alloc + D_EXTPAD;
d = x_alloc + D_EXTPAD + half + D_EXTPAD + D_EXTPAD;
memcpy(r, x, half * sizeof(float));
memcpy(d, x + half, half * sizeof(float));
for (i = 1; i <= D_EXTPAD; i++)
{
r[-i] = r[i];
r[(half - 1) + i] = r[half - i];
d[-i] = d[i - 1];
d[(half - 1) + i] = d[(half - 1) - i];
}
for (n = half; n--;)
{
*x++ = h_b[0] * r[0] + h_b[2] * (r[1] + r[-1]) + g_b[3] * (d[1] + d[-2]) +
g_b[1] * (d[0] + d[-1]);
*x++ = h_b[1] * (r[1] + r[0]) + h_b[3] * (r[2] + r[-1]) + g_b[4] * (d[2] + d[-2]) +
g_b[2] * (d[1] + d[-1]) + g_b[0] * d[0];
d++;
r++;
}
}
void f97_2D(float **rows, short width, short height, short levels, short inverse)
{
short x, y, w, h, l;
long yw;
float *buffer, *rows_ptr;
/*Check the dimensions for compatability. */
if (width % (1 << levels) || height % (1 << levels))
{
printf("width and height must be divisible by 2^levels");
exit(10);
}
if ((x_alloc =malloc(sizeof(float) * (width + height + F_EXTPAD + F_EXTPAD))) == NULL)
{
printf("malloc failed");
exit(10);
}
/* Allocate a work array (for transposing columns) */
if ((buffer = malloc(sizeof(float) * height)) == NULL)
{
printf("malloc failed");
exit(10);
}
/*Compute the rows wavelet transform. */
if (!inverse)
{ /*forward transform. */
for (l = 0; l < levels; l++)
{
w = width >> l;
h = height >> l;
/*Rows. */
for (y = 0; y < h; y++)
forwardf97(rows[y], w);
/* Columns. */
rows_ptr = rows[0];
for (x = 0; x < w; x++)
{
for (y = 0, yw = 0; y < h; y++, yw += width)
buffer[y] = rows_ptr[yw];
forwardf97(buffer, h);
for (y = 0, yw = 0; y < h; y++, yw += width)
rows_ptr[yw] = buffer[y];
rows_ptr++;
}
}
} else
{
for (l = levels - 1; l >= 0; l--)
{
w = width >> l;
h = height >> l;
/* Columns. */
rows_ptr = rows[0];
for (x = 0; x < w; x++)
{
for (y = 0, yw = 0; y < h; y++, yw += width)
buffer[y] = rows_ptr[yw];
inversef97(buffer,
