commit : denoise object added

This commit is contained in:
Rafat Hussain 2017-09-24 18:24:45 +05:30
parent a26e271971
commit f906ae8d66
7 changed files with 290 additions and 44 deletions

View File

@ -1,6 +1,26 @@
#include "denoise.h" #include "denoise.h"
denoise_object denoise_init(int length, int J,char* wname) {
denoise_object obj = NULL;
obj = (denoise_object)malloc(sizeof(struct denoise_set) +sizeof(double));
obj->N = length;
obj->J = J;
strcpy(obj->wname,wname);
//Set Default Values
strcpy(obj->dmethod,"sureshrink");
strcpy(obj->ext,"sym");
strcpy(obj->level,"first");
strcpy(obj->thresh,"soft");
strcpy(obj->wmethod,"dwt");
return obj;
}
void visushrink(double *signal,int N,int J,char *wname,char *method,char *ext,char *thresh,char *level,double *denoised) { void visushrink(double *signal,int N,int J,char *wname,char *method,char *ext,char *thresh,char *level,double *denoised) {
int filt_len,iter,i,dlen,dwt_len,sgn, MaxIter,it; int filt_len,iter,i,dlen,dwt_len,sgn, MaxIter,it;
double sigma,td,tmp; double sigma,td,tmp;
@ -25,8 +45,6 @@ void visushrink(double *signal,int N,int J,char *wname,char *method,char *ext,ch
dwt(wt,signal); dwt(wt,signal);
} else if(!strcmp(method,"swt")) { } else if(!strcmp(method,"swt")) {
swt(wt,signal); swt(wt,signal);
} else if(!strcmp(method,"modwt")) {
modwt(wt,signal);
} else { } else {
printf("Acceptable WT methods are - dwt,swt and modwt\n"); printf("Acceptable WT methods are - dwt,swt and modwt\n");
exit(-1); exit(-1);
@ -102,8 +120,6 @@ void visushrink(double *signal,int N,int J,char *wname,char *method,char *ext,ch
idwt(wt,denoised); idwt(wt,denoised);
} else if(!strcmp(method,"swt")) { } else if(!strcmp(method,"swt")) {
iswt(wt,denoised); iswt(wt,denoised);
} else if(!strcmp(method,"modwt")) {
imodwt(wt,denoised);
} }
free(dout); free(dout);
@ -259,3 +275,85 @@ void sureshrink(double *signal,int N,int J,char *wname,char *method,char *ext,ch
wave_free(wave); wave_free(wave);
wt_free(wt); wt_free(wt);
} }
void denoise(denoise_object obj, double *signal,double *denoised) {
if(!strcmp(obj->dmethod,"sureshrink")) {
sureshrink(signal,obj->N,obj->J,obj->wname,obj->wmethod,obj->ext,obj->thresh,obj->level,denoised);
} else if(!strcmp(obj->dmethod,"visushrink")) {
visushrink(signal,obj->N,obj->J,obj->wname,obj->wmethod,obj->ext,obj->thresh,obj->level,denoised);;
} else {
printf("Acceptable Denoising methods are - sureshrink and visushrink\n");
exit(-1);
}
}
void setDenoiseMethod(denoise_object obj, char *dmethod) {
if (!strcmp(dmethod, "sureshrink")) {
strcpy(obj->dmethod, "sureshrink");
}
else if (!strcmp(dmethod, "visushrink")) {
strcpy(obj->dmethod, "visushrink");
}
else {
printf("Acceptable Denoising methods are - sureshrink and visushrink\n");
exit(-1);
}
}
void setDenoiseWTMethod(denoise_object obj, char *wmethod) {
if (!strcmp(wmethod, "dwt")) {
strcpy(obj->wmethod, "dwt");
}
else if (!strcmp(wmethod, "swt")) {
strcpy(obj->wmethod, "swt");
}
else {
printf("Wavelet decomposition method can be either dwt or swt");
exit(-1);
}
}
void setDenoiseWTExtension(denoise_object obj, char *extension) {
if (!strcmp(extension, "sym")) {
strcpy(obj->ext, "sym");
}
else if (!strcmp(extension, "per")) {
strcpy(obj->ext, "per");
}
else {
printf("Signal extension can be either per or sym");
exit(-1);
}
}
void setDenoiseParameters(denoise_object obj, char *thresh,char *level) {
//Set thresholding
if (!strcmp(thresh, "soft")) {
strcpy(obj->thresh, "soft");
}
else if (!strcmp(thresh, "hard")) {
strcpy(obj->thresh, "hard");
}
else {
printf("Thresholding Method - soft or hard");
exit(-1);
}
// Set Noise estimation at the first level or at all levels
if (!strcmp(level, "first")) {
strcpy(obj->level, "first");
}
else if (!strcmp(level, "all")) {
strcpy(obj->level, "all");
}
else {
printf("Noise Estimation at level - first or all");
exit(-1);
}
}
void denoise_free(denoise_object object) {
free(object);
}

View File

@ -11,11 +11,38 @@ Copyright (c) 2017, Rafat Hussain
extern "C" { extern "C" {
#endif #endif
//depends on J typedef struct denoise_set* denoise_object;
denoise_object denoise_init(int length, int J,char* wname);
struct denoise_set{
int N; //signal length
int J; // Levels of Wavelet decomposition
char wname[10]; //Wavelet name
char wmethod[10]; //Wavelet decomposition method - dwt or swt
char ext[10]; // Signal Extension - sym or per
char thresh[10]; // thresholding - soft or hard
char level[10]; // Noise Estimation level - first or all
char dmethod[20]; //Denoising Method -sureshrink or visushrink
//double params[0];
};
void visushrink(double *signal,int N,int J,char *wname,char *method,char *ext,char *thresh,char *level,double *denoised); void visushrink(double *signal,int N,int J,char *wname,char *method,char *ext,char *thresh,char *level,double *denoised);
void sureshrink(double *signal,int N,int J,char *wname,char *method,char *ext,char *thresh,char *level,double *denoised); void sureshrink(double *signal,int N,int J,char *wname,char *method,char *ext,char *thresh,char *level,double *denoised);
void denoise(denoise_object obj, double *signal,double *denoised);
void setDenoiseMethod(denoise_object obj, char *dmethod);
void setDenoiseWTMethod(denoise_object obj, char *wmethod);
void setDenoiseWTExtension(denoise_object obj, char *extension);
void setDenoiseParameters(denoise_object obj, char *thresh,char *level);
void denoise_free(denoise_object object);
#ifdef __cplusplus #ifdef __cplusplus
} }

View File

@ -11,6 +11,35 @@ int compare_double(const void* a, const void* b)
} }
double mean(double* vec, int N) {
int i;
double m;
m = 0.0;
for (i = 0; i < N; ++i) {
m+= vec[i];
}
m = m / N;
return m;
}
double var(double* vec, int N) {
double v,temp,m;
int i;
v = 0.0;
m = mean(vec,N);
for (i = 0; i < N; ++i) {
temp = vec[i] - m;
v+= temp*temp;
}
v = v / N;
return v;
}
double median(double *x, int N) { double median(double *x, int N) {
double sigma; double sigma;
@ -242,3 +271,52 @@ void getDWTRecCoeff(double *coeff,int *length,char *ctype,char *ext, int level,
free(out); free(out);
} }
void autocovar(double* vec,int N, double* acov,int M) {
double m,temp1,temp2;
int i,t;
m = mean(vec,N);
if ( M > N) {
M = N-1;
printf("\n Lag is greater than the length N of the input vector. It is automatically set to length N - 1.\n");
printf("\n The Output Vector only contains N calculated values.");
} else if ( M < 0) {
M = 0;
}
for(i = 0; i < M;i++) {
acov[i] = 0.0;
for (t = 0; t < N-i;t++) {
temp1 = vec[t] - m;
temp2 = vec[t+i] - m;
acov[i]+= temp1*temp2;
}
acov[i] = acov[i] / N;
}
}
void autocorr(double* vec,int N,double* acorr, int M) {
double var;
int i;
if (M > N) {
M = N - 1;
printf("\n Lag is greater than the length N of the input vector. It is automatically set to length N - 1.\n");
printf("\n The Output Vector only contains N calculated values.");
}
else if (M < 0) {
M = 0;
}
autocovar(vec,N,acorr,M);
var = acorr[0];
acorr[0] = 1.0;
for(i = 1; i < M; i++) {
acorr[i] = acorr[i]/var;
}
}

View File

@ -19,9 +19,12 @@
extern "C" { extern "C" {
#endif #endif
int compare_double(const void* a, const void* b); int compare_double(const void* a, const void* b);
double mean(double* vec, int N);
double var(double* vec, int N);
double median(double *x, int N); double median(double *x, int N);
double mad(double *x, int N); double mad(double *x, int N);
@ -35,6 +38,10 @@ void getDWTDetail(wt_object wt, double *detail, int N, int level);
void getDWTRecCoeff(double *coeff,int *length,char *ctype,char *ext, int level, int J,double *lpr, void getDWTRecCoeff(double *coeff,int *length,char *ctype,char *ext, int level, int J,double *lpr,
double *hpr,int lf,int siglength,double *reccoeff); double *hpr,int lf,int siglength,double *reccoeff);
void autocovar(double* vec,int N,double* acov, int M);
void autocorr(double* vec,int N,double* acorr, int M);
#ifdef __cplusplus #ifdef __cplusplus
} }

View File

@ -1,30 +0,0 @@
/*
Copyright (c) 2017, Rafat Hussain
*/
#ifndef DENOISE_H_
#define DENOISE_H_
#include "wavelib.h"
#ifdef __cplusplus
extern "C" {
#endif
//depends on J
void visushrink(double *signal,int N,int J,char *wname,char *method,char *ext,char *thresh,char *level,double *denoised);
void sureshrink(double *signal,int N,int J,char *wname,char *method,char *ext,char *thresh,char *level,double *denoised);
void getDWTRecCoeff(double *coeff,int *length,char *ctype,char *ext, int level, int J,double *lpr,
double *hpr,int lf,int siglength,double *reccoeff);
double mad(double *x, int N);
#ifdef __cplusplus
}
#endif
#endif /* DENOISE_H_ */

56
header/wauxlib.h Normal file
View File

@ -0,0 +1,56 @@
/*
Copyright (c) 2017, Rafat Hussain
*/
#ifndef WAUXLIB_H_
#define WAUXLIB_H_
#include "wavelib.h"
#ifdef __cplusplus
extern "C" {
#endif
typedef struct denoise_set* denoise_object;
denoise_object denoise_init(int length, int J,char* wname);
struct denoise_set{
int N; //signal length
int J; // Levels of Wavelet decomposition
char wname[10]; //Wavelet name
char wmethod[10]; //Wavelet decomposition method - dwt or swt
char ext[10]; // Signal Extension - sym or per
char thresh[10]; // thresholding - soft or hard
char level[10]; // Noise Estimation level - first or all
char dmethod[20]; //Denoising Method -sureshrink or visushrink
//double params[0];
};
void visushrink(double *signal,int N,int J,char *wname,char *method,char *ext,char *thresh,char *level,double *denoised);
void sureshrink(double *signal,int N,int J,char *wname,char *method,char *ext,char *thresh,char *level,double *denoised);
void denoise(denoise_object obj, double *signal,double *denoised);
void setDenoiseMethod(denoise_object obj, char *dmethod);
void setDenoiseWTMethod(denoise_object obj, char *wmethod);
void setDenoiseWTExtension(denoise_object obj, char *extension);
void setDenoiseParameters(denoise_object obj, char *thresh,char *level);
void denoise_free(denoise_object object);
void getDWTRecCoeff(double *coeff,int *length,char *ctype,char *ext, int level, int J,double *lpr,
double *hpr,int lf,int siglength,double *reccoeff);
double mad(double *x, int N);
#ifdef __cplusplus
}
#endif
#endif /* WAUXLIB_H_ */

View File

@ -2,7 +2,8 @@
#include <stdlib.h> #include <stdlib.h>
#include <string.h> #include <string.h>
#include <math.h> #include <math.h>
#include "../header/denoise.h"
#include "../header/wauxlib.h"
static double rmse(int N,double *x,double *y) { static double rmse(int N,double *x,double *y) {
double rms; double rms;
@ -49,7 +50,9 @@ int main() {
// gcc -Wall -I../header -L../Bin denoisetest.c -o denoise -lwauxlib -lwavelib -lm // gcc -Wall -I../header -L../Bin denoisetest.c -o denoise -lwauxlib -lwavelib -lm
double *sig,*inp,*oup; double *sig,*inp,*oup;
int i,N,J; int i,N,J;
FILE *ifp,*ofp; FILE *ifp;
denoise_object obj;
double temp[2400]; double temp[2400];
char *wname = "db5"; char *wname = "db5";
@ -100,16 +103,23 @@ int main() {
for(i = 0; i < N;++i) { for(i = 0; i < N;++i) {
inp[i] = temp[i]; inp[i] = temp[i];
} }
visushrink(inp,N,J,wname,method,ext,thresh,level,oup); obj = denoise_init(N,J,wname);
setDenoiseMethod(obj,"visushrink");// sureshrink is also the default. The other option is visushrink
setDenoiseWTMethod(obj,method);// Default is dwt. the other option is swt
setDenoiseWTExtension(obj,ext);// Default is sym. the other option is per
setDenoiseParameters(obj,thresh,level);// Default for thresh is soft. Other option is hard
// Default for level is first. The other option is all
denoise(obj,inp,oup);
// Alternative to denoise_object
// Just use visushrink and sureshrink functions
//visushrink(inp,N,J,wname,method,ext,thresh,level,oup);
//sureshrink(inp,N,J,wname,method,ext,thresh,level,oup); //sureshrink(inp,N,J,wname,method,ext,thresh,level,oup);
//ofp = fopen("denoiseds.txt", "w"); //ofp = fopen("denoiseds.txt", "w");
for(i = 0; i < N;++i) {
//fprintf(ofp,"%g \n",oup[i]);
}
//fclose(ofp);
printf("RMSE %g\n",rmse(N,sig,inp)); printf("RMSE %g\n",rmse(N,sig,inp));
printf("Corr Coeff %g\n",corrcoef(N,sig,inp)); printf("Corr Coeff %g\n",corrcoef(N,sig,inp));
@ -119,6 +129,6 @@ int main() {
free(sig); free(sig);
free(inp); free(inp);
free(oup); denoise_free(obj);
return 0; return 0;
} }