#include <stdio.h>
#include <stdlib.h>
int powerof2(int exponent) {
if ( exponent < 0 ) {
return 0;
}
if ( exponent == 0 ) {
return 1;
}
return 2 << (exponent - 1);
}
/*---------------------------------------------------------------------------*/
/*
Calculate one iteration of the Daubechies forward FWT in 1-dimension.
*/
void Daubechies_forward_pass_1d (int n, float * s)
{
int i;
int npts;
float * a = NULL;
float * c = NULL;
const float h[4] = { 0.683013, 1.18301, 0.316987, -0.183013 };
npts = powerof2 (n);
a = (float *) malloc (sizeof(float) * npts/2);
c = (float *) malloc (sizeof(float) * npts/2);
for (i = 0; i < npts/2; i++)
{
a[i] = (h[0]*s[(2*i)%npts] + h[1]*s[(2*i+1)%npts] + h[2]*s[(2*i+2)%npts]
+ h[3]*s[(2*i+3)%npts]) / 2.0;
c[i] = (h[3]*s[(2*i)%npts] - h[2]*s[(2*i+1)%npts] + h[1]*s[(2*i+2)%npts]
- h[0]*s[(2*i+3)%npts]) / 2.0;
}
for (i = 0; i < npts/2; i++)
{
s[i] = a[i];
s[i + npts/2] = c[i];
}
free (a); a = NULL;
free (c); c = NULL;
}
/*---------------------------------------------------------------------------*/
/*
Calculate the Daubechies forward fast wavelet transform in 1-dimension.
*/
void Daubechies_forward_FWT_1d (int n, float * s)
{
int m;
int npts;
npts = powerof2 (n);
for (m = n-1; m >= 0; m--)
{
Daubechies_forward_pass_1d (m+1, s);
printf(npts,s[m]); //for (i=0;i<32;i++) printf("x[%d]=%f\n",i,x[i]);
printf("\n");
}
}
int main()
{
float x[]={1,2,3,4,2,23,12,3,2};
Daubechies_forward_FWT_1d(2,x);
return 0;
}
