# 一，一維信號的拓延
在Mallat算法的推導中,假定輸入序列是無限長的,而實際應用中常常是分時采樣,即輸入序列為有限長.此時,濾波器系數與輸入序列卷積時就會出現輪空的現象.因此有必要對原始信號進行邊界延拓,減小邊界誤差.解決的方法通常有補零法和周期延拓法.
1)補零法是在輸入序列的末尾補零.補零法的缺點是可能人為造成輸入序列邊界的不連續,從而使得較高頻率的小波系數很大.
2)周期延拓法是將原來有限長的輸入序列拓展成周期的序列.周期延拓可適用于任何小波變換,但可能導致輸入序列邊緣的不連續,使得高頻系數較大.這種方式的拓延卷積后與源信號的長度一致。
3)對稱延拓（matlab默認采取這種方式）可避免輸入序列邊緣的不連續,但只適用于對稱小波變換.本文根據Mallat算法公式,編寫了對稱延拓方式的小波變換的一般實現方法.
注：筆者采用的編譯器為VS2013，當前系統為win10。
### 1，常見的3種邊緣拓延方法
設輸入信號為f(n),長度為srcLen,濾波器長度為filterLen.下面給出信號邊界處理幾種方法的具體表達式如下：
1）周期拓延：
? ? ? ? ? ??f(n+ srcLen) ? , ?-(filterLen -1)≤ n< 0
? f′(n)=??f(n) ? ? ? ? ? ? ? ? , ?0≤ n≤ srcLen?-1
? ? ? ? ? ? f(n -srcLen) ? ?, ?srcLen -1< n≤srcLen+filterLen?-2
舉例說明：以“1 2 3 4 5 6 7 8”這個長度為8的信號為例，當濾波器的長度為4時，其具體的拓延長度為6（單邊為3）：
6 7 8 （1 2 3 4 5 6 7 8）1 2 3
2）對稱延拓（本文重點）
? ? ? ? ? ? f(-n -1) ? ? ? ? ? ? ? ? ? ? , ? -(filterLen-1)≤ n< 0
f′(n)=? ? f(n) ? ? ? ? ? ? ? ? ? ? ? ? ?, ? 0≤ n≤srcLen-1
? ? ? ? ? ? f(2srcLen -n -1) ? ? ? , ? srcLen-1< n≤ srcLen+ filterLen-2
舉例說明：以“1 2 3 4 5 6 7 8”這個長度為8的信號為例，當濾波器的長度為4時，其具體的拓延長度為6（單邊為3）：
3 2 1 （1 2 3 4 5 6 7 8）8 7 6
3）零值填補
? ? ? ? ? ? ?0 ? ? ? ? ? ? ? ?, ? -(filterLen -1)≤ n< 0
f′(n)=?? ? f(n) ? ? ? ? ? ?, ?0≤ n≤srcLen?-1
? ? ? ? ? ? ?0 ? ? ? ? ? ? ? ?, ?srcLen?-1< n≤srcLen+filterLen?-2
舉例說明：以“1 2 3 4 5 6 7 8”這個長度為8的信號為例，當濾波器的長度為4時，其具體的拓延長度為6（單邊為3）：
0 0 0 （1 2 3 4 5 6 7 8）0 0 0?
#二，Mallat算法分解過程
在Mallat算法中,信號分解過程按照系數形式,
其中:hd(n)和gd(n)分別為分解的低通和高通濾波器系數,長度為filterLen?.
吐槽：信號濾波的過程其實就是濾波器與信號卷積的過程，也就是濾波器的頻譜和信號的頻譜相乘的過程。
那么低通濾波器頻譜是高頻低，低頻高，與信號相乘時高頻的會被過濾掉；高通濾波器同理！
### 1，詳細的分解過程
以db2小波為例，通過對稱拓延后的詳細計算過程如下：
從上圖可知，小波的Mallat算法分解后的序列長度由原序列長srcLen和濾波器長filterLen決定。從Mallat算法（采用對稱拓延）的分解原理可知，分解后的序列就是原序列與濾波器序列的卷積再進行隔點抽取而來。即分解抽取的結果長度為(srcLen+filterLen-1)/2。
上述方法也是matlab默認采取的方法，在MATLAB中輸入代碼，查看結果便知：
~~~
[Lo_D,Hi_D,Lo_R,Hi_R] = wfilters('db2');%獲取小波系數
xx=[420.2, 423.53, 423.52, 423.35, 424.52, 428, 430.79, 428.92];
wname='db2';
level=1;
[c,l]=wavedec(xx,level,wname);%進行一層分解
~~~
結果為（被存儲在系數矩陣c中，l記錄的是每段系數的長度）：
595.429871699852??597.374655846484598.449909371632?604.108998389031607.245626013478
-2.03920021086643 ??0.509501871423169-1.289053089583672.33989974581161??-1.14513645475084
與上圖一致
### 2，C++實現分解算法
獲取db2小波系數：
~~~
//系數精度很高，采自matlab
extern double db2_Lo_D[4] = { -0.129409522550921, 0.224143868041857, 0.836516303737469, 0.482962913144690 };
extern double db2_Hi_D[4] = { -0.482962913144690, 0.836516303737469, -0.224143868041857, -0.129409522550921 };
extern double db2_Lo_R[4] = { 0.482962913144690, 0.836516303737469, 0.224143868041857, -0.129409522550921 };
extern double db2_Hi_R[4] = { -0.129409522550921, -0.224143868041857, 0.836516303737469, -0.482962913144690 };
~~~
C++實現如下：
~~~
// 一維信號的小波分解
int CWavelet::DWT(
double *pSrcData,//分解的源信號
int srcLen,//源信號的長度
double *pDstCeof//分解出來的,本函數將返回此長度
)
{
//本程序禁止出現這種情況，否則數據出錯（對稱拓延長度為filterLen-1，如果大于了signalLen將越界）
if (srcLen < m_dbFilter.filterLen - 1)
{ //實際上信號的長度可以是任意的（matlab順序：信號拓延-》卷積-》下采樣），
//但是本程序為了算法速度，寫法上不允許
cerr << "錯誤信息：濾波器長度大于信號!" << endl;
Sleep(1000);
exit(1);
}
int exLen = (srcLen + m_dbFilter.filterLen - 1) / 2;//對稱拓延后系數的長度
int k = 0;
double tmp = 0.0;
for (int i = 0; i < exLen; i++)
{
pDstCeof[i] = 0.0;
pDstCeof[i + exLen] = 0.0;
for (int j = 0; j < m_dbFilter.filterLen; j++)
{
k = 2 * i - j + 1;
//信號邊沿對稱延拓
if ((k<0) && (k >= -m_dbFilter.filterLen + 1))//左邊沿拓延
tmp = pSrcData[-k - 1];
else if ((k >= 0) && (k <= srcLen - 1))//保持不變
tmp = pSrcData[k];
else if ((k>srcLen - 1) && (k <= (srcLen + m_dbFilter.filterLen - 2)))//右邊沿拓延
tmp = pSrcData[2 * srcLen - k - 1];
else
tmp = 0.0;
pDstCeof[i] += m_dbFilter.Lo_D[j] * tmp;
pDstCeof[i + exLen] += m_dbFilter.Hi_D[j] * tmp;
}
}
return 2 * exLen;
}
~~~
處理結果為：（與matlab一致）
附帶測試主函數：
~~~
int main()
{
system("color 0A");
double s[8] = { 420.2, 423.53, 423.52, 423.35, 424.52, 428, 430.79, 428.92};
int signalLen = sizeof(s) / sizeof(double);
cout << "原始信號：" << endl;
for (int i = 0; i < signalLen; i++)
cout << s[i] << " ";
cout << endl;
CWavelet cw(2);
double *dst = new double[10];
cw.DWT(s, signalLen,dst);
cout << "分解后的系數：" << endl;
for (int i = 0; i < 10; i++)
cout << dst[i] << " ";
cout << endl;
delete[] dst;
dst = NULL;
system("pause");
return 0;
}
~~~
# 三，Mallat算法重構過程
### 1，詳細的重構過程
以db2小波為例，通過對稱拓延后的詳細重構過程如下：
在matlab中我們可以通過以下代碼實現重構（完全實現原信號）
~~~
[Lo_D,Hi_D,Lo_R,Hi_R] = wfilters('db2');%獲取小波系數
xx=[420.2, 423.53, 423.52, 423.35, 424.52, 428, 430.79, 428.92];
wname='db2';
level=1;
[c,l]=wavedec(xx,level,wname);%進行一層分解
cc=waverec(c,l,wname);
~~~
###2，C++實現重構算法
~~~
// 一維小波反變換，重構出源信號
void CWavelet::IDWT(
double *pSrcCoef,//源分解系數
int dstLen,//重構出來的系數的長度
double *pDstData//重構出來的系數
)
{
int p = 0;
int caLen = (dstLen + m_dbFilter.filterLen - 1) / 2;
for (int i = 0; i < dstLen; i++)
{
pDstData[i] = 0.0;
for (int j = 0; j < caLen; j++)
{
p = i - 2 * j + m_dbFilter.filterLen - 2;
//信號重構
if ((p >= 0) && (p<m_dbFilter.filterLen))
pDstData[i] += m_dbFilter.Lo_R[p] * pSrcCoef[j] + m_dbFilter.Hi_R[p] * pSrcCoef[j + caLen];
}
}
}
~~~
處理結果（與matlab一致）
附帶測試主函數：
~~~
int main()
{
system("color 0A");
double s[8] = { 420.2, 423.53, 423.52, 423.35, 424.52, 428, 430.79, 428.92};
int signalLen = sizeof(s) / sizeof(double);
cout << "原始信號：" << endl;
for (int i = 0; i < signalLen; i++)
cout << s[i] << " ";
cout << endl;
CWavelet cw(2);
double *dst = new double[10];
cw.DWT(s, signalLen,dst);
cout << "分解后的系數：" << endl;
for (int i = 0; i < 10; i++)
cout << dst[i] << " ";
cout << endl;
double *dstsrc = new double[signalLen];
cw.IDWT(dst, signalLen, dstsrc);
cout << "重構后的系數：" << endl;
for (int i = 0; i < signalLen; i++)
cout << dstsrc[i] << " ";
cout << endl;
delete[] dstsrc;
dstsrc = NULL;
delete[] dst;
dst = NULL;
system("pause");
return 0;
}
~~~
#參考資源：
【1】喬世杰.小波圖像編碼中的對稱邊界延拓法[ J].中國圖像圖形學報,2000,5(2):725-729.
