贴一个小波变换算法的C源码,大家指点下!
#define N0 128#include "stdio.h"
#include "stdlib.h"
#include "math.h"
#include "string.h"
void db4(double *h,double *g,double *hh,double *gg);
void wd(int N,double *h,double *g,double *c0,double *c,double *d);
void wr(int N,double *h,double *g,double *c, double *d,double *cd);
void main()
{
double fk,c0,c,d;
double h,g,hh,gg;
float fk0;
FILE *fp;
int i,k,j,n,l,N;
fp=fopen("wdata.dat","rt");
fscanf(fp,"%d",&N);
for(k=0;k<N;k++) fscanf(fp,"%f",&fk0);
fclose(fp);
db4(h,g,hh,gg);
for(k=0;k<N;k++) {
c0=fk0;
c=0;
d=0;
}
wd(N,hh,gg,c0,c,d);
wr(N,hh,gg,c,d,c0);
for(k=0;k<N;k++) printf("k=%d c0=%f c=%f\n",k,fk0,c0);
return;
}
void wd(int N,double *h,double *g,double *c0,double *c,double *d)
/* wavelet decomposition */
{
int k,n,k2,l;
double ck,dk;
for(k=0;k<N;k++) {
ck=0.0;
dk=0.0;
for(l=0;l<8;l++) {
n=k+l;
ck+=c0*h;
dk+=c0*g;
}
c=ck;
d=dk;
}
for(k=0;k<N/2;k++) {
k2=2*k;
c0=c;
c0=d;
}
return;
}
void wr(int N,double *h,double *g,double *c,double *d,double *c0)
/* wavelet reconstruction */
{
int k,n,l,k2;
double ck,cn,dn;
for(k=0;k<N/2;k++) {
k2=2*k;
c=c0;
c=0;
d=c0;
d=0;
}
for(k=0;k<N;k++) c0=0.0;
for(k=0;k<N;k++) {
ck=0.0;
for(l=0;l<8;l++) {
n=k-l;
cn=c[(N+n)%N];
dn=d[(N+n)%N];
ck+=cn*h+dn*g;
}
c0=ck;
}
return;
}
void db4(double *h,double *g,double *hh,double *gg)
/* Daubechies 4 wavelet */
{
int k,isgn;
h=-0.0105974017850890;
h= 0.0328830116668852;
h= 0.0308413818355607;
h=-0.1870348117190931;
h=-0.0279837694168599;
h= 0.6308807679398597;
h= 0.7148465705529154;
h= 0.2303778133088964;
isgn=1;
for(k=0;k<8;k++) {
gg=isgn*h;
isgn=-isgn;
}
for(k=0;k<8;k++) {
g=gg;
hh=h;
}
return;
}
float fun(float x)
{
float pi=3.1415926;
float yx=30*exp(-x/40)*sin(2*pi*x/40);
return(yx);
}
//这个用的是 DB4小波，周期延拓，可以实现精确重构的
[ 本帖最后由 风花雪月 于 2008-6-17 15:01 编辑 ] 我回去验证一下
不知道您的计算周期是多少
我编的程序计算周期不满足要求
验证完毕后回来再顶你的贴
[ 本帖最后由 wr_9864 于 2008-3-22 23:44 编辑 ]
程序可读性差！
现在很多个人写的程序都是这样，缺乏必要的注释过一段时间可能自己也忘记了 没看懂。。。。。。。。。。。。。。。。。。
回复 楼主 的帖子
一维的，而且没有采用提升算法不知道我理解的对不？？ 没有看懂呀，如果加上注释就好了
页:
[1]