一、多分辨率分析的图像处理
如果物体的尺寸较小或对比度较低,通常以较高的分辨率来研究;
如果物体的尺寸较大或对比度较高,粗略的观察就已足够
如果较大物体和较小物体同时存在时,以不同的分辨率来研究它们将更具优势,即多分辨率分析。
1、图像金字塔
图像金字塔是一系列以金字塔形状排列的,分辨率逐步降低的图像集合,底部是一个高分辨率图像,顶部是低分辨率图像。
金字塔低分辨率级别用于分析较大的结构或图像的整体内容;高分辨率图像适合于分析单个物体的特性。
预测残差:一幅图像与该图像的一个预测版本间的差,预测残差可更有效率地编码。
预测残差金字塔:
步骤1:第j级图像通过近似滤波器,进行步长为2的下采样(近似滤波器若是邻域平均,会产生平均金字塔;若是高斯滤波,会产生高斯金字塔;若是不进行过滤,会产生子取样金字塔)
步骤2:把步骤1的结果通过插值(最近邻,双线性,双三次内插)进行上采样
步骤3:把输入图像与步骤2得到的图像做差,形成预测残差
2、子带编码
一幅图像被分解为一组频带受限的分量称为子带。
子带可以重组在一起无失真地重建原始图象
每个子带通过对输入进行带通滤波而得到
子带带宽小于原始图像带宽,子带可以进行无信息损失的抽样
原始图象的重建可以通过内插、滤波、和叠加单个子带来完成
分析滤波器组即为分解,综合滤波器组即为重构。
分解就是通过下采样提取图像的低频(背景)和高频(细节)信息。每一层(尺度)的分解都是对上一层分解中的低频信息进行再分解。
重构是通过上采样将分解的低频和高频信息再合并为分解前的图像的近似。
3、哈尔变换
哈尔变换的基函数是一直最古老最简单的正交小波。
二、多分辨率展开
1、级数展开
一个信号或函数f(x)可展开为函数的线性组合:
k为有限和或无限和的整数下标,αk为展开系数,φk(x)为展开函数。
若展开唯一,即任意给定的f(x)只有一组αk与之对应。则φk(x)称为基函数,集合{φk(x)}称为可表示这样一类函数的基。
可展开的函数形成了一个函数空间,称为展开集合的闭合跨度,表示如下:
对任意函数空间V及相应的展开集合{φk(x)},都有一个表示为{φk ̃(x)}的对偶函数集合。展开系数αk即是对偶函数φk ̃(x)和函数f(x)的内积。
2、尺度函数
考虑由实、平方可积函数φ(x)的整数平移和二值尺度组成的展开函数集合{φj,k(x)}
其中φj,k(x)对所有的j,k∈Z(整数集)和φ(x)属于L2®都成立(L2®表示度量的、平方可积的一维函数集合,R为实数集)。
整数平移k决定了φj,k(x)沿x轴的位置;尺度j决定了φj,k(x)的宽度;2j/2控制函数的幅度。
因为φj,k(x)的形状随j发生变化,故φ(x)称为尺度函数。
选择适当的φ(x),可使{φj,k(x)}张成L2®。
对任意的j,将k上张成的子空间表示为:
增加j就会增加Vj的大小,进而允许子空间中包含具有更小变量或更细细节的函数。因为随着j的增大,用于表示子空间函数的φj,k(x)会变窄,且x有较小变化就可以分开。
hφ(n)为小波函数系数,hφ为小波向量。
三、小波变换
1、小波级数展开
与小波ψ(x)和尺度函数φ(x)相关的函数f(x)∈L2®的小波级数展开如下:
j0是任意的开始尺度,cj0(k)和dj(k)分别是尺度函数和小波函数下的展开系数αk的改写形式。cj0(k)称为近似和/或尺度系数,dj(k)称为细节和/小波系数。
即为被展开的函数和展开函数的内积。
2、离散小波变换
上一小节的小波系数展开将一个连续变量函数映射为一系列系数。若待展开的函数是离散的(即数字序列),则得到的系数就称为离散小波变换(DWT)。则序列f(n)的正向DWT系数(类似于上一节的cj0(k)和dj(k))如下:
其中φj0,k(n)ψj,k(n)是基函数φj0,k(x)ψj,k(x)的取样形式,n = 0,1,2,…,M-1。所以f(n)的展开如下:
3、连续小波变换
离散小波变换的自然延伸是连续小波变换(CWT),连续小波变换将一个连续函数变换为两个连续变量(平移和尺度)的高冗余度函数。得到的变换易于解释并且对于时间-频率分析是有价值的。我们感兴趣的是离散图像,此处只是为了完整性。
连续平方可积函数f(x)的连续小波变换与实数值小波ψ(x)的关系定义如下:
其中
s和τ分别为尺度参数和平移参数。
f(x)可由连续小波反变换求得:
其中
式中Ψ(μ)是ψ(x)的傅里叶变换。
四、快速小波变换
快速小波变换(FWT)是一种实现离散小波变换(DWT)的高效计算。类似于2子带的子带编码方案。
正变换,即分解
先定义尺度函数和小波函数,再计算尺度j的近似和细节系数。
3.2.1 一维情况
由前面描述的尺度函数和小波函数
其中hφ(n)为尺度函数系数,hφ为尺度向量;hψ(n)为小波函数系数,hψ为小波向量。
用2j对x进行尺度化,用k对它进行平移,并令m = 2k + n,则
尺度函数:
小波函数
近似系数和细节系数如下:
hφ(-n) = h0(n)(2带子带编码和解码系统中的分析部分),hψ(-n) = h1(n)
反变换,即重构
快速小波反变换
反变换即综合滤波器,综合滤波器和分析滤波器之间的顺序是相反的。
g0(n) =h0(-n) = hφ(n) =,g1(n) = h1(-n) = hψ(n)
五、 二维小波变换
在二维情况下,需要一个二维尺度函数φ(x,y)和三个二维小波ψH(x,y),ψV(x,y),ψD(x,y)。每个二维小波都是两个一维函数的乘积。排除产生一维结果的乘积,如φ(x)ψ(x),4个剩下的乘积可产生尺度函数
和可分的“方向敏感”小波
这些小波度量函数的变化——图像的灰度变化——沿不同方向的变化:ψH度量沿列方向的变化(如水平边缘),ψV响应沿行方向的变化(如垂直边缘),ψD度量对应对角线方向的变化。
定义尺度和平移基函数:
故大小为M*N的图像f(x,y)的离散小波变换是:
通常令j0 = 0,并选择N = M = 2J,故j = 0,1,2,…,J-1和m = n = 0,1,2,…,2j-1。
f(x,y)可通过离散小波反变换得到:
利用小波去噪:
import matplotlib.pyplot as plt
import pywt
# Get data:
ecg = pywt.data.ecg() # 生成心电信号
index = []
data = []
for i in range(len(ecg)-1):
X = float(i)
Y = float(ecg[i])
index.append(X)
data.append(Y)
# Create wavelet object and define parameters
w = pywt.Wavelet('db8') # 选用Daubechies8小波
maxlev = pywt.dwt_max_level(len(data), w.dec_len)
print("maximum level is " + str(maxlev))
threshold = 0.04 # Threshold for filtering
# Decompose into wavelet components, to the level selected:
coeffs = pywt.wavedec(data, 'db8', level=maxlev) # 将信号进行小波分解
for i in range(1, len(coeffs)):
coeffs[i] = pywt.threshold(coeffs[i], threshold*max(coeffs[i])) # 将噪声滤波
datarec = pywt.waverec(coeffs, 'db8') # 将信号进行小波重构
mintime = 0
maxtime = mintime + len(data) + 1
plt.figure()
plt.subplot(2, 1, 1)
plt.plot(index[mintime:maxtime], data[mintime:maxtime])
plt.xlabel('time (s)')
plt.ylabel('microvolts (uV)')
plt.title("Raw signal")
plt.subplot(2, 1, 2)
plt.plot(index[mintime:maxtime], datarec[mintime:maxtime-1])
plt.xlabel('time (s)')
plt.ylabel('microvolts (uV)')
plt.title("De-noised signal using wavelet techniques")
plt.tight_layout()
plt.show()