#include <stdio.h>
#include <stdlib.h>
#define LENGTH 512//信号長度
/******************************************************************
* 一維卷積函數
*
* 說明: 循環卷積,卷積結果的長度與輸入信号的長度相同
*
* 輸入參數: data[],輸入信号; core[],卷積核; cov[],卷積結果;
* n,輸入信号長度; m,卷積核長度.
*
* 李承宇, [email protected]
*
* 2010-08-18
******************************************************************/
void Covlution(double data[], double core[], double cov[], int n, int m)
{
int i = 0;
int j = 0;
int k = 0;
//将cov[]清零
for(i = 0; i < n; i++)
{
cov[i] = 0;
}
//前m/2+1行
i = 0;
for(j = 0; j < m/2; j++, i++)
{
for(k = m/2-j; k < m; k++ )
{
cov[i] += data[k-(m/2-j)] * core[k];//k針對core[k]
}
for(k = n-m/2+j; k < n; k++ )
{
cov[i] += data[k] * core[k-(n-m/2+j)];//k針對data[k]
}
}
//中間的n-m行
for( i = m/2; i <= (n-m)+m/2; i++)
{
for( j = 0; j < m; j++)
{
cov[i] += data[i-m/2+j] * core[j];
}
}
//最後m/2-1行
i = (n - m) + m/2 + 1;
for(j = 1; j < m/2; j++, i++)
{
for(k = 0; k < j; k++)
{
cov[i] += data[k] * core[m-j-k];//k針對data[k]
}
for(k = 0; k < m-j; k++)
{
cov[i] += core[k] * data[n-(m-j)+k];//k針對core[k]
}
}
}
/******************************************************************
* 一維小波變換函數
*
* 說明: 一維小波變換,隻變換一次
*
* 輸入參數: input[],輸入信号; output[],小波變換結果,包括尺度系數和
* 小波系數兩部分; temp[],存放中間結果;h[],Daubechies小波基低通濾波器系數;
* g[],Daubechies小波基高通濾波器系數;n,輸入信号長度; m,Daubechies小波基緊支集長度.
*
* 李承宇, [email protected]
*
* 2010-08-19
******************************************************************/
void DWT1D(double input[], double output[], double temp[], double h[],
double g[], int n, int m)
{
// double temp[LENGTH] = {0};//?????????????
int i = 0;
/*
//尺度系數和小波系數放在一起
Covlution(input, h, temp, n, m);
for(i = 0; i < n; i += 2)
{
output[i] = temp[i];
}
Covlution(input, g, temp, n, m);
for(i = 1; i < n; i += 2)
{
output[i] = temp[i];
}
*/
//尺度系數和小波系數分開
Covlution(input, h, temp, n, m);
for(i = 0; i < n; i += 2)
{
output[i/2] = temp[i];//尺度系數
}
Covlution(input, g, temp, n, m);
for(i = 1; i < n; i += 2)
{
output[n/2+i/2] = temp[i];//小波系數
}
}
void main()
{
double data[LENGTH];//輸入信号
double temp[LENGTH];//中間結果
double data_output[LENGTH];//一維小波變換後的結果
int n = 0;//輸入信号長度
int m = 6;//Daubechies正交小波基長度
int i = 0;
char s[32];//從txt檔案中讀取一行資料
static double h[] = {.332670552950, .806891509311, .459877502118, -.135011020010,
-.085441273882, .035226291882};
static double g[] = {.035226291882, .085441273882, -.135011020010, -.459877502118,
.806891509311, -.332670552950};
//讀取輸入信号
FILE *fp;
fp=fopen("data.txt","r");
if(fp==NULL) //如果讀取失敗
{
printf("錯誤!找不到要讀取的檔案/"data.txt/"/n");
exit(1);//中止程式
}
while( fgets(s, 32, fp) != NULL )//讀取長度n要設定得長一點,要保證讀到回車符,這樣指針才會定位到下一行?回車符傳回的是零值?是,非數字字元經過atoi變換都應該傳回零值
{
// fscanf(fp,"%d", &data[count]);//一定要有"&"啊!!!最後讀了個回車符!适應能力不如atoi啊
data[n] = atof(s);
n++;
}
//一維小波變換
DWT1D(data, data_output, temp, h, g, n, m);
//一維小波變換後的結果寫入txt檔案
fp=fopen("data_output.txt","w");
//列印一維小波變換後的結果
for(i = 0; i < n; i++)
{
printf("%f/n", data_output[i]);
fprintf(fp,"%f/n", data_output[i]);
}
//關閉檔案
fclose(fp);
}
轉載于:https://www.cnblogs.com/hjj-fighting/p/11106862.html