2009-7-19 1
第六章光 学 小 波 变 换第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
2
第六章 光学小波变换
6,1 引言
6,2 从短时傅里叶变换到小波变换
6,3 小波变换的定义和性质
6,4 实现一维小波变换的光学系统
6,5 用多通道匹配滤波实现二维小波变换
6,6 光学小波变换匹配滤波器在图像识别中的应用
6,7 光学 Haar 小波变换和图形边缘探测
6.1 引 言如果 g(x)是一个时域或空域中分布在 (-?,?)中的稳恒过程或稳定分布,则傅里叶分析给出了近乎完美的结果,然而,在自然界和科学技术中还有大量信号,它们具有局部的或定域的特性,例如语言信号,声纳信号,
各种电脉冲等 。 这些 信号只出现在一个短暂的时间问隔内,此后很快衰减到零,称快速过程或暂态过程,一个很短暂的信号,可以 称为,小波,信号对于局部信号或暂态过程,傅里叶分析就不完全适用,首先,我们仅对?t内的时间信号感兴趣,没有必要在过去,现在及未来的无限长时间范围内对信号进行分析；类似地,在处理定域于?x 内的空间图像时,
也没有必要对全平面内的信号进行全面的分析,
在许多情况下,在?t 或?x 以外的信号是未知的,
它可能是零,也可能是背景噪声；对它们我们不太了解,测不准,或不感兴趣,如不加选择地把 (-?,?)内全部信号进行傅里叶处理,还可能产生较大的误差甚至错误,此外,一个局部的信号在?t 或?x 以外较远处几乎完全等于零,当用它们的频谱来恢复或重构这些信号时,在?t 或?x 外很远处也会出现一些非零的分量,它们一般不是信号,而是 在傅里叶逆变换中频域综合不够充分而产生的噪声,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
5
在一些课题中,我们往往不满足于了解信号在全部区间内的综合的频谱分布,而 希望了解某一区间或某些区间内信号对应的频谱,例如在地震勘探中,为了分辨分层的地层和矿床结构,我们需要在时域和频域中仔细分析不同时刻的信号在不同频谱区间中的行为,而傅里叶分析只能提供在长时间内的信号整体的频谱,显然不能满足我们的要求,
近年来发展起来的小波分析,正好克服了傅里叶分析的上述缺点,它和傅里叶分析的一个重要区别,在于它恰恰适用于处理局部或暂态信号,因此,小波分析成为信号分析,图像处理,
数据压缩,语音信号分析等领域中的重要工具；
在地震勘探信号处理,边缘探测,语音信号合成中则有特殊的用途,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
6
6.2 从短时傅里叶变换到小波变换
6.2.1 短时傅里叶变换 (STFT)
为了有效地提取一个局部信号 g(x)的信息,引入一个局部化的变换,
所谓局部化,包含两个要素：
第一,被分析的区间要有一定的宽度
x,我们仅对?x及其附近的信息进行处理；
第二,被分析的区间有一个中心坐标 xc，
当 xc 改变时,就可以提取不同的信息,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
7
为了实现局部化,一个有效的方案是在傅里叶变换中加一个 窗函数 w(x)：
Gw(v,xo) =?∞- ∞g(x)exp(-i2?vx) w(x- xo) dx (1)
= [W(v)exp(-i2?vxo)]* G(v) (2)
式中 W 和 G 分别是 w 和 g 的傅里叶变换,
只要 w(x)和 W(v)有足够快的衰减速度,窗函数就是一个局部化的函数,
函数 f (x) 和 g(x) 的 内积的定义,
( f (x),g(x) ) =?∞-∞ f (x) g(x)dx
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
8
窗函数的中心 xc定义,
xc = (w(x),x w(x))/(w(x),w(x))
=?∞- ∞ w(x) xw(x)dx /?∞- ∞ w(x)w(x)dx
式中 (*,*)表示函数 f 和 g 的内积,xc与 xo不一定相等。
窗函数的宽度 则定义：
w = 2[(w(x),(x- xc )2 w(x))/(w(x),w(x))]1/2
= 2[?∞-∞w(x)(x-xc )2w(x)dx /?∞-∞ w(x) w(x)dx]1/2
注意,它是节 1.2.4中所定义的信号空域宽度的两倍 。
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
9
Gw(v,xo) =?∞- ∞g(x)exp(-i2?vx) w(x- xo) dx (1)
由于窗函数具有局部处理的功能,因此 (1)式定义的变换称为 短时傅里叶变换 (STFT)．
STFT和 FT的一个重要区别,
频率变量 v 和坐标变量 xo同时出现在变换函数中,
在 STFT中,窗口宽度则隐含于 Gw(v,xo) 内,
正是 xo和窗口宽度?w,使这一变换具有局部处理的功能,改变 xo,窗口就在空域中移动,以获取不同区域的信息,xo 通常称为 位移因子 ；?w 则限制了被处理空间的范围,
频率窗中心,vc = (W(v),vW(v))/(W(v),W(v)) (5)
频率窗宽度,
W=2[(W(v),(v- vc)2W(v)w(x))/(W(v),W(v))]1/2 (6)
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
10
当?w和?W 都有限时，我们称函数 w(x)在空域和频域同时局部化．
w?W称为 空间 -频率窗,它限制了空域和频域中被处理区域的范围。根据?w,?W的定义及测不准关系式，并注意信号宽度定义的区别，我们就有,?w?W≥1/? (7)
当 高斯函数取为窗函数 时,(7)式中的等式成立,这种情况下 STFT具有最小处理区域,
STFT的局部性,其特征在于处理过程限制在空间 -频率窗内进行,且窗的位置是可变的,然而无论?w还是?W都是常数,不会随信号中心频率的变化而变化,这使 STFT在处理一些奇异性的信号时显得无力,恰恰是在这一点上,小波变换 具备比 STFT更强的功能,
6.2.2 Gabor 变换早在 1946年,Gabor就提出了下面的变换
(8)
称 Gabor变换,其中?和 b为变换的参数,上式又可表为
(9)
式中
(10)
因此 Gabor变换就是高斯窗短时傅里叶变换,
窗函数中心坐标 xc = 0,窗的宽度?w = 1.414?
w?(x)的傅里叶变换,W?(v) = exp(-2?2?2v2) (13)
也是高斯函数,频率窗宽度?W = 1/ 1.414 (14)
因此有?w?W = 1 /?,
图中将空域和频域同时表达出来,称 空间 -频率坐标系,空 -频窗 则表示为图中的一个矩形,Gabor变换空 -频窗的高度和宽度都是恒定的 。
Gabor变换在频域中的表达式,[式中? = 1/(2)2 ]
可见 Gabor变换在频域和空域中的表达式具有相似的形式．
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
13
Gabor变换的特点,
(1)实现空域和频域处理的局部化 [中心位于 (b,v)
空 -频窗为 1.414?× (1/ 1.414) = 1/? ]
(2) 由 (8)式和 (15)式,
可看出变换是参数?,b 和变量 v 的函数。上式给出的 积分是一个调制包络,载波 exp(-i2?vb)的频率 (即中心频率 ) v 与参数?无关,不会随? 的变化而变化，这正是所有短时傅里叶变换共同的 缺点 。
(17)
6.2.3 Morlet小波变换为了克服 Gabor变换中窗口尺寸不能变动的缺点,
Gabor变换的 基元函数 可改写为子函数定义信号函数 g(x)的 Morlet小波变换,
称 变换的母函数,
引入参数 a,b，
Morlet小波变换与 Gabor变换的实质性差别,
小波变换,vc= v/a,?w =1.414?a,?W=1/1.414a,
当 vc 增高时 (a 减小 ),?w 变小而?W增大，可处理更多的高频信息；当 vc 降低时 (a 增大 ),?W变小而
w加宽，可容纳足够多个空间周期，以保证处理精度．
Morlet小波变换,在处理低频信号时空间窗自动加宽,
在空间窗范围内包含的信号空间周期相同,这就保证了小波变换以同样的精度去处理不同中心频率的信号,这正是小波变换与短时傅里叶变换的根本区别,
Gabor变换,窗的宽度是常数，当 vc增高时，一定宽度的空间窗内包含的空间周期增加，所以 变换的精度是随频率而变化的 ；
小波变换的空间-
频率窗第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
16
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
17
6.3 小波变换的定义和性质
6.3.1 小波变换的定义母函数 h(x)的基本小波函数 ha,b(x)定义为
(1)
式中 b称为小波变换的 位移因子,a > 0 称为 伸缩因子,上式表明基本小波是母函数经平移和缩放的结果,基本小波又简称 小波,
信号函数 g(x)的 小波变换定义 为,(2)
由于相关运算较易用光学相关器进行,因此小波变换可以用我们已熟悉的光学相关系统来实现,
)a bx(ha1)x(h b,a
)b(g)ab(ha1dx)x(g)a bx(ha1)x(g),x(h)}x(g{W *b,ab,a
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
18
Morlet小波的 母函数 是子函数 是
Mor1et小波的基函数式中,m,n = 0,1,2,… ； a = aom ； b = nbo.
)2 te x p ()t2ie x p (2 1)t(h 2
2
0
)a bt(ha1)t(h n,m
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
19
当 m = n = 0时,a = 1,b = 0,即为母函数；
当 m = 1,n = 0时,a = a0,b = 0,对应的一阶小波函数为当 m = -1,n = 0时,a = a0-1,b = 0,有负一阶小波函数
h(t),h1,0(t),和 h-1,0(t)已分别在上图中画出,
)a2 te x p ()a t2ie x p (
2
1
a
1)t(h
2
0
2
2
0
0
0
0,1
)2 tae x p ()ta2ie x p (
2
1
a/1
1)t(h
2
22
0
00
0
0,1
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
20
并非任何函数都可以作为小波变换的函数 h(x),h(x)必须在?x时衰减到零,实际使用的小波变换母函数 h(x),当
x 时迅速衰减,使它的不显著为零的分量只存在于一个很小的区间内,这正是,小波,名称的来由 。 实际上,也只有迅速衰减的小波才使变换 (2)式具备局部化的特征,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
21
6.3.2 小波变换在频域中的表达式在频域中,小波可表示为式中 H(v)是小波母函数 h(x)的傅里叶变换,在空域的扩大 x/a 等价于频域的压缩 av,空域的位移 b 等价于频域的位相移动 exp(-i2?vb)．
由 Parseval定理,小波变换的傅里叶变换为
(4)
其中 H 和 G 分别是 h 和 g 的傅里叶变换,上式表明 信号 g(x)的小波变换可以用 4f系统实现,
)a(H)b2ie x p (adx)x(h)x2ie x p ()x(H b,ab,a
d)(G)b2ie x p ()a(Ha)}x(g{W *b,a?
6.3.3 逆变换和相容性条件小波变换 (2)式的逆变换定义为其中 Ch满足条件
(6)式称,相容性条件,,是逆变换存在的条件,
下面我们来证明它,ha,b(x) 可表为以 (4),(7)式代入 (5)式,得到上式中最内部的积分为?(v’ - v),因此上式成立的条件是即 (6)式．当 v = 0时,相容性条件 要求,H(0) = 0,(10)
即小波函数没有零频分量．由于 H(0)=?-h(x)dx = 0 (11)
意味着 h(x)必须是振荡函数,平均值为零,其傅里叶谱直流分量为零,
6.3.4 正则性从理论上讲,任何满足相容性条件的函数都可当作小波变换的母函数,然而在实用中,为了使变换具备局部化的功能,h(x)和 H(v)在空域和额域中都是迅速衰减的,
它们不显著为零的分量分别分布于空域和领域中的原点附近,
此外，要求 Wa,b作为 a的函数，应当是充分光滑的，
当 a?0 时,Wa,b?0，即要求 Wa,b在 a=0附近是正则的,设
b = 0，则有将 g(x)在 x =0的邻域内展开成泰勒级数：
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
25
代入 (12)式得到式中此是小波函数 h 的 n 阶矩：
Mn =?- h(?)?nd? (n = 0,1,…) (15)
由 (11)式 Mo = H(0) = 0 (16)
设 Mp = 0 ( p = 0,1,…,n) (17)
则在 0 的邻域内第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
26
随着 a? 0,Wa.o? 0 的速率为即对于一个足够平滑的函数 g(x),Wa.o{g}以 a n+1/2
的速率随 a 趋近于零，称它为 n 阶小波函数,
由节 1.1(24)式,Mp = 0 ( p = 0,1,…,n)意味着 H(p) (v) = 0 ( p = 0,1,…,n)
表明 v = 0是 H(v) 的 n 阶零点,
此外,(16)式 (即相容性条件 )保证 Wa.o随 a 趋于零的速度的下限为 a1/2．
6.3.5 小波变换的空间 -频率窗和处理过程的局部化小波变换在空域中的处理局限于空间窗内
(23)
小波变换在频域中的处理局限于频率窗内
(29)
空间 -频率的处理就局限于空间 -频率窗内：
(31)
小波变换处理过程的特点：
(1)空间窗宽度 a?w 和频率窗宽度?W/a 均随 a 的变化而变化,窗的面积?w?W 与 a 无关,
(2) 中心频率 v /a 与带宽 (即频率窗宽 )之比
Q = (v /a)/ (?W/a) = v /?W (32)
与中心频率大小无关,仅取决于 H(v),Q是测量精度的特征量,上式表明 小波变换的测量精度与频率无关,
当 v/a增大时 (a 减小 ),频率窗自动变宽,使小波变换在不同频率下具有相同的检测精度；反之,当 v/a 减小时 (a 增大 ),空间窗自动加宽,以容纳同样数目的信号空间周期,有人把小波变换的这种性能比喻为,自动变焦,(zooming)，
伸缩因子 a 常称小波变换的 频率变量,位移因子 b 则称为 坐标变量,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
29
6.3.6 常用小波函数
1,Haar小波
Haar小波是双极性阶跃函数它以 t = 1/2为中心的奇对称实函数,满足小波存在条件,
)]43t(2[r e c t)]41t(2[r e c t)t(h
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
30
Haar小波的傅里叶变换 是它的模是正的偶函数，以 v = 0 为对称轴,相位因子 exp(-i?v)是由于 h(t)是以 t =1/2为中心的奇对称性所引起的．
Haar小波对于离散的缩放因子和位移量是正交的，其傅里叶谱的振幅 |H(v)|随 1/v 很慢地收敛到零．
/)]c os (1)[ie xp (i2)(H vvvv
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
31
2,Morlet小波
Morlet小波是由分析声像技术引人的,
可表示为其实部是余弦 -高斯函数,
Morlet小波的傅里叶谱 是平移到 vo 及 –vo
处的两个高斯函数，即显然,H(v)是正的实偶函数．
)2te x p ()t2ie x p ()t(h
2
0
]})(2e x p [])(2{ e xp [2)(H 202202
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
32
图,h(t)的实部和 h(t) 的傅里叶谱由于 H(0)≧ 0,所以小波变换的存在条件没有满足,但是对于较大的 vo,H(0)十分接近于 0。 在数值计算时可近似地看作为 0.
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
33
3,Mexican – hat 小波这种型如墨西哥帽的小波函数被广泛用于零交叉多分辨边缘检测,其母函数实际上是高斯函数的二阶导数,即该函数是实偶函数，满足小波变换存在条件,
Mexican - hat小波的傅里叶变换是它也是实偶函数,
)
2
t
e x p ()t1()t(h
2
2
)2e x p (4)(H 222
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
34
Mexican - hat小波函数及其傅里叶变换高斯函数的高阶导数也可以用来作为小波函数,高斯函数 n 阶导数的傅里叶谱 H(v)
将是乘以 ( i 2?v )n 的高斯函数,所以 H(0) =
0，满足子波变换的存在条件。 Mexican –
hat 小波变换收敛速度很快．
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
35
4,Meyer 小波
Meyer小波是用其傅里叶变换来定义的,即
H(v) = exp(-i?v )sin[ω(v)]
式中,ω(v)是偶对称函数，如图所示．
图中,AB弧段有一个对称中心，在 v =1/2 和
ω(1/2) =?/4 处,即
3231),(2)1(
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
36
而 BC弧段具有相同的形状，是 AB弧段的翻转和展开，其对称中心在 v = 1 和 ω(1) = π/4，即如果 H(v)中相移因子 ex(-i?v)被忽略，则 H(v)为实偶函数．相移因子的作用，是使 h(t)在时间轴上的位移为 t =1/2．
Meyer小波可表示 为它也是一个以 t =1/2 为对称点的实函数，衰减较快,Meyer已经证明，具有离散缩放和位移因子的这种小波，可以构成一组正交基底．
3231),(2)2(
d])21t(2c o s [)](s i n [2)t(h
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
37
6.3.7 离散小波变换 (DWT)
正如傅里叶变换一样,小波变换也可以分为两种形式：连续的和离散的变换,
对于 离散小波变换,我们将假定信号 g(x)也是离散的,离散小波变换的形式为对应 离散小波逆变换 为
a/]a/)xx[(h)x(gW
x
jjiij
a/]a/)xx[(hW)x(s
i j
jjiij
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
38
6.4 实现一维小波变换的光学系统从小波变换的定义可知,N维信号函数的小波变换是 2N维函数,因此计算的工作量很大,尽管专门用于小波变换的超大规模集成电路已有报道,但人们仍在考虑用光学系统来实现小波变换,
因为光学信息处理器具有高度的并行处理性能,
一维小波变换光学处理系统第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
39
xy-L2-uv 构成 x 方向的一维傅里叶变换器。在
SLM1上输入信号 g(x)，经 L2 的变换作用，在 uv
平面上形成它的傅里叶谱 G(u)．
在 uv 平面上放置第二个空间光调制器 SLM2，
它被分成 M个沿 u方向的带状区域，这些带状区域中分别显示具有不同伸缩因子 am 的基元函数 h
的傅里叶谱：
H*(amu) (m = 1,2,…,M),(1)
假定 H 是实的,从而有
H*(amu) = H(amu),(2)
{H(amu)| m = 1,2,…,M } 构成多通道小波变换匹配滤波器,G(u) 经滤波后成为
H*(amu) G(u) (m = 1,2,…,M),(3)
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
40
xy-L2-uv 构成 x 方向的一维傅里叶变换器。在
SLM1上输入信号 g(x)，经 L2 的变换作用，在 uv
平面上形成它的傅里叶谱 G(u)．
在 uv 平面上放置第二个空间光调制器 SLM2，
它被分成 M个沿 u方向的带状区域，这些带状区域中分别显示具有不同伸缩因子 am 的基元函数 h
的傅里叶谱：
H*(amu) (m = 1,2,…,M),(1)
假定 H 是实的,从而有
H*(amu) = H(amu),(2)
{ H(amu)| m = 1,2,…,M } 构成多通道小波变换匹配滤波器,G(u) 经滤波后成为
H*(amu) G(u) (m = 1,2,…,M),(3)
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
41
uv - L3 -构成像散系统在 子午面 (vz平面 )内,L3使 uv 平面成像在上 ;
在 弧矢面 (uz平面 )内,柱面镜没有作用,uz 位于球面镜的前焦面, 位于其后焦面,构成沿 u 方向的一维傅里叶逆变换,参见图 6.8．
图 6.8(a)在子午面内构成成像系统
(b)在弧矢面内构成一维傅里叶逆变换系统第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
42
由于成像作用,在平面上沿?方向相应形成 uv 平面上各带状通道的像,对于第 m 个通道,
由于沿 u 方向的傅里叶逆变换作用,得到
(a,am) =? H*(amu) G(u) exp(i2?a u) (4)
由节 6.3(4)式,有
(5)
在图像处理系统中,将 CCD输出的信号除以
(am )1/2即得到小波变换,
在 中,伸缩因子 a 是分立的,由一组 M个滤波器 引入处理器；而位移因子
b 则是连续的,与输入平面的坐标 a成正比,
) M,2,1,m ( )a,a(a}g{W mma,a m
}g{W a,a m
),a(H m u
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
43
光学 Morlet 小波变换实例实验采用下图所示的光学多通道小波变换相关器．两个柱面透镜的焦距为 240 mm，球面透镜焦距为 210 mm．相关器由 He-Ne激光作光源．
实现一维小波变换的二维光学相关器小波函数是余弦 ——高斯 Morlet 函数它的傅里叶变换是
H(v)是一个实函数，表达式
H(0) = exp(-2?2?2v02)≠ 0，
但是随余弦函数中的频率 v0 的增加而迅速下降．
例如,当 2?v0 = 5,σ= 1时,H(0) = 3.7× 10-6，
因此小波变换的存在条件在一个非常好的近似范围内被满足．
2 te x p)t2c o s (
2
1)t(h
2
2
0
]})(2e x p [])(2{ e x p [21)(H 20222022
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
45
光学多通道小波变换相关器中光学小波变换的实验结果第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
46
6.5 用多通道匹配滤波实现二维小波变换
6.5.1 单通道小波变换系统一维小波变换定义,
二维小波变换也可类似定义：
(1)
为了简单起见，令 ax = ay = a，即 x,y方向按相同的尺度加以缩放，得到
(2)
)b(g)ab(ha1dx)x(g)a bx(ha1)}x(g{W *b,a
dd),(
a
b
,
a
b
h
aa
1
)b,b(
a
b
,
a
b
h
aa
1
)}y,x({W
y
y
x
x*
yx
yx
y
y
x
x
yx
b,b,a,a yxyx
dd),(a b,a bha1 )}y,x({W yx*b,b,a yx
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
47
在频域中,(1)式变成 匹配滤波的频域表达式 (3)
匹配滤波可以用标准的 4f系统实现,如图 6.9所示,
图 6.9 用 4f 系统实现小波坐换将二维信号函数?(x,y)经过 SLMl 输入系统,
则在 Ll的频谱面上将出现它的谱?(u,v)．
d)dbb(2ie x p [),()a,a(Haa}{W yxyx*yxb,b,a,a yxyx vuvuvuvu
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
48
在谱面上放置第二个 SLM2,将匹配滤波函数 H*(axu,ayv)通过 SLM2 对?(u,v) 进行滤波,
则形成 H*(axu,ayv)?(u,v),再经过第二个透镜
L2,在输出平面上得到它的傅里叶逆变换,
由上式可知,H*?的傅里叶逆变换,即信号?(x,
y) 的小波变换
d)dbb(2ie x p [),()a,a(Haa}{W yxyx*yxb,b,a,a yxyx vuvuvuvu
)}y,x({W yxyx b,b,a,a?
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
49
我们注意到,位移因子 ( bx,by )是与输出平面的坐标对应的变量,但伸缩因子
( ax,ay )却是给定的,亦即我们只能对给定的伸缩因子 ( ax,ay ) 实现小波变换,不同的 ( ax,ay )的变换只能通过依次输入不同的匹配滤波函数 H*( axu,ayv ) 来实现,速度很慢,发挥不了光学系统并行处理的优越性,
6.5.2 用 Dammann 光栅进行多通道相关处理由于 comb(x)与任意函数?(x)相乘,等于对?(x) 的抽样； comb(x)与?(x)的乘积再进行傅里叶变换,在谱面上得到 comb(u)与?(u)的卷积,相当于谱项?(u)在 comb(u)
中各个? 函数的位置上重复出现,只要?(x) 的带宽是有限的,则我们总可以通过足够密集的抽样手续,使各谱项在频域内互相分离,即 对于给定的输入信号?(x,y),可利用梳状函数在谱面上复制出一系列谱函数?mn(u,v)，
m,n = 0,± 1,± 2,…,如果每个?mn 代表一个处理通道,就可以实现多通道的并行处理,在这里,关键的器件是梳状函数器件,又称 Dammann光栅,
一维 Dammann光栅的透过率函数可表为
(4)
式中?是栅线宽度,d 是空间周期,见图 6.10(a)．
mdxr e c tdxc o m b*xr e c td1)x(t
m
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
51
若将 Dammann光栅放在 4f系统的输入面上,则在谱面上得到
T(u) =? sinc(?u)comb(du),(5)
它是一个被 sinc函数调制的梳状函数，见图
6.0(b),若?足够小，则 sinc(?u) 变化缓慢，
在中心附近的一些?函数幅度接近于 1．当然,? 越小，透过的光能量就越少．
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
52
把两个同样规格的一维 Dammann光栅相对旋转 90o并叠在一起,就得到二维 Dammann光栅,其透过率函数可表为
(6)
在?很小时,上式可近似表为
(7)
显见 二维 Dammann光栅是一个? 函数点阵,它们分布在间隔为 d 的栅线的交点上,这些地方的透过率为 1,其余地方的透过率为 0．
dyc o m b*yr e c tdxc o m b*xr e c td 1)y,x(t 2?
dy,dxc o m bd1)y,x(t 2?
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
53
在上图所示的 4f系统的输入平面上放置信号
(x,y),把二维 Dammann光栅紧贴在?(x,y) 上，
在光栅后的光场复振幅分布为
(x,y) comb(x/d,y/d) /d2 (8)
在 uv 平面上得到它的傅里叶谱
(9)
式中 fo = 1/d,(10)
d
1)nf,mf(
d
1
)dd(c o m b*)()(
m n
n,m2
m n
oo2
vu
vu,vu,vu,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
54
m.n =?(u - mfo,v - nfo) (11)
是一系列中心位于 Amn=(mfo,nfo) 处复现的谱项,如 设计一个匹配滤波器列阵,位于 Amn
处的复数透过率为 H*m.nexp[i2?(upm+vqn)]，
亦即：
(12)
滤波器中各谱项是按角度编码的，相当于用不同方向传播的平面波 exp [ i2?(upm+vqn)]
作为参考光制成的全息匹配滤波器，通常是计算机生成的全息匹配滤波器 (CGH)．
e)]nf(a),mf(a[H)(F
m n
)qp(2
omom
* nm vuvuvu,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
55
经滤波后,再通过 4f 系统中第二个傅里叶透镜 L2的作用,在输出平面上得到傅里叶逆变换：
式中 Cmn 是常数系数,(13)
m n
qb,pb,a,amn2
m n
nymx
n
ny
m
mx
mn2
m n
)]qb()pb([2i
nm
*
mn2
)]qb()pb([2i
oo
m n
onom
*
2yx
)}y,x({WC
d
1
)qb,pb(
a
qb
,
a
pb
hC
d
1
dde),()a,a(HC
d
1
dde)nf,mf(
)]nf(a),mf(a[H
d
1
)b,b(
nymxnm
nymx
nymx
vuvuvu
vuvu
vu
vu
vu
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
56
亦即频域中的相移形成输出平面上的位移，以 ( pm,qn )为中心形成一系列在空间相互分离的项，每一项都代表一个不同的伸缩因子 ( am,an )的小波变换，位移因子则由以 ( pm,qn )为中心的坐标来表示,
这样一来，我们就用 多通道 4f 系统实现了二维光学小波变换，它是缩放因子的分立函数，是位移因子的连续函数．
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
57
6.5.3用体全息存储器进行多通道相关处理重铬酸明胶 (DCG) 或 光 折 变 晶 体
(BaTiO3,LiNbO3,BSO等 )都可以做成体全息存储器,并可以用角度对全息图进行编码,这一性质也可以用于光学小波变换,
参见图 6,11．
将缩放后的母函数 hmn = h( x/am,y/an ) (14)
经 SLM1输入系统,则在频谱面上产生了它的傅里叶变换
Hmn = am an H ( amu,anv ),在谱平面上放置体全息记录器件,例如光折变晶体,并用倾斜的平行参考光
Rmn = exp[ -i2?(upm+vqn)] (15)
照射,与物光 am an H ( amu,anv ) 相干叠加,形成全息图,一系列 { hmn } 依次输入系统,并用不同角度的参考光 { Rmn } 编码,最后得到,(16)
式中 第三项正是我们所需的小波变换匹配滤波函数,
eHaaeHaaHaa1{
RHaa)(F
m n
)qp(2
mnnm
)qp(2*
mnnm
2*
mn
2
n
2
m
m n
2
mn
*
mnnm
nmnm
vuvu
vu,
当我们将所分析的信号?(x,y)通过 SLM1输入系统,
它的谱?(u,v) 经过 F(u,v)的滤波后得到
(17)
我们略去了与匹配滤波无关的项,再经 L2 的傅里叶逆变换,在输出平面上得到同样得到以 (pm,qn)为中心的一系列空间互相分离的小波变换．
)e()aa(Haa
m n
)qp(2
nm
*
nm
nm vuvu,v u,
m n
qb,pb,a,a
m n
nymx
n
ny
m
mx
)]qb()pb([2i
m n
nm
*
nmyx
)}y,x({W
)qb,pb(
a
qb
,
a
pb
h
dde),()a,a(Haa)b,b(
nymxnm
nymx
vuvuvu
vu
(18)
(16)式前二项经逆变换重合于原点，但第四项则是以 (-pm,-qn) 为中心的卷积项，在输出平面上与相关项关于原点对称，可称为,鬼像,
(ghost image)，因此在角度编码时要避免这些鬼像对其他小波变换项的干扰．
目前全息存储密度已做得相当高，因此在晶体中可以存储相当多个按方向编码的小波变换匹配滤波器．
eHaaeHaaHaa1{
RHaa)(F
m n
)qp(2
mnnm
)qp(2*
mnnm
2*
mn
2
n
2
m
m n
2
mn
*
mnnm
nmnm
vuvu
vu,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
61
6.6 光学小波变换匹配滤波器在图像识别中的应用图像识别或特征识别,
从大量信息或背景中检测某一特定的图形或指定的特征信息,并排斥其他图形信息 。
在一般情况下,我们只知道需要检测的图形的特征,对于其他图形的特征,我们事先可能并不知道,或知之甚少,但检测系统必须排斥这些图形,图像识别系统的这种性能称为,排他性,。
光学小波变换识别系统在这方面的性能比常规的相关识别系统更强,因此有可能利用这一效应设计成有应用价值的小波变换图像识别系统 。
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
62
6.6.1 边缘增强效应图形的重要特征之一是它的 形状或轮廓 。 为了识别某一特定的图形,往往只需认定它的轮廓,而并不需要研究它的内部细节 。 轮廓就是图形的边缘,一旦图形的边缘被清晰地勾画出来,这一图形就容易识别了 。 相对于图形整体而言,边缘显然是局部,因此我们可以期望小波变换在边缘检测中有特殊的功效 。
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
63
选择,墨西哥帽,式母函数,定义为
(1)
引入 g(x,y) = exp[-(x2 + y2)/2] (2)
则有 h(x,y) = -(1/2?)?2g(x,y) (3)
由 F {?2g(x,y)} = - 4?2(u2+ v2)G(u,v),得到
H(u,v) = - 2? (u2+ v2)G(u,v) (4)
式中 G 是 g 傅里叶变换
G(u,v) = 2? exp[ - 2? (u2+ v2)] (5)
代入 (4)得到：
H(u,v) = 4?2(u2+ v2) exp[ - 2? (u2+ v2)] (6)
2
yxe x p)]yx(2[
2
1)y,x(h 2222
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
64
图 6.12给出 a = 0.6 及 a = 0.3 的墨西哥帽小波函数及它们的傅里叶谱,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
65
函数?(x,y) 的小波变换
7)：
亦即? 的小波变换是? 和缩放后的高斯函数 g 的相关的二次导数,
),(
a
y
,
a
x
g
a2
1
-
dd),(
a
y
,
a
x
g
a2
1
-
dd),(
a
y
,
a
x
h
a
1
)y,x(
a
y
,
a
x
h
a
1
}{W
2
2
*
y,x,a
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
66
由于
(8)
伸缩因子 a 正是高斯函数的特征尺度 。 g 和
相关的结果则是平滑效应,? 中比 a 小得多的精细结构都被平滑掉 。
此外,? 和 g 的相关运算的结果再二次求导,?(x,y) 中振幅不变的区域 (常数项 )
及线性变化的区域 (一次项 )都等于零,而在振幅变化的拐点两侧不为 0。 边界正是这样的拐点,
a2
yxe x p
a
y,
a
xg
2
22
),(ay,axga2 1-}{W 2y,x,a?
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
67
(a) 边界函数?(x)和高斯函数 g(x)
(b)Cr (a,x)= g(x/a)(x) (c)小波变换函数 D (x)
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
68
图 (a)：边界函数?(x),高斯函数 g(x/a)；
图 (b)：给出它们的相关,相关运算对?(x)
起到了平滑的作用,结果噪声都被平均掉；
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
69
可以 在边界内外侧看到小波变换的一对正、负峰，
它们明确指示了边界的位置,边界内、外区域内小波变换函数都是零．其最终效果恰恰是边界突出或轮廓突出，这正是我们所期望的．
在二维小波变换的情况下，伸缩因子 a是矢量而不是标量,a = (am,an)，一般情况下 am ≠an,
从而在两个方向上的平滑效果可分别加以控制．
图 (c),给出小波变换作为位移因子的函数 (差一个常数因子 )，即
2[g(x,y)(x)
待处理坦克模型图像和小波函数平面图像及小波变换边缘增强后实验结果图 3,待处理图像小波函数图像及小波变换边缘增强后图像
( d ) I = 2 4,a = 0,0 8
( e ) I = 2 4,a = 0,1 5
( c ) I = 2 4,a = 0,0 5
(f) I = 2 4,a = 0,2 0
( a ) 待处理图像
( b ) 小波函数图像待检测彩色图像灰度图像后边缘检测
RG
B
空间边缘检测
YUV
空间边缘检测待检测彩色图像灰度图像后边缘检测
Lαβ
空间边缘检测选择域值二值化处理第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
73
6.6.2 小波变换匹配滤波图像识别课题,一般情况下仅仅勾画出图像或区域的轮廓或边界是不够的,还要求认定输入图像中是否包含要求识别的目标,
在常规光学信息处理中,我们用匹配滤波方法达到这一目的,在小波变换信息处理中,匹配滤波方案同样适用,
设 输入 要求识别的 目标图像?(x,y) 和 另一图像? (x,y),首先对它进行二维小波变换，可采用上一节中介绍的各种方法，其结果，得到了边缘增强的小波变换谱函数 Wa,x,y{?} 和?的小波变换
Wa,x,y{?}，然后再以 W{?}和 W{?} 作为输入信号,进入第二个光学相关识别系统,例如 4f 系统，
或第五章中介绍的各种实现傅里叶变换的系统。
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
74
在谱面上,W{?} 和 W{?} 的频谱函数分别由以下二式给出：
(10)
(11)
表明：在第二个相关识别系统中,当输入函数为
W{?} 和 W{?} 时,谱面上将出现?H* 和?H*。
我们 用?(x,y)所对应的小波变换匹配滤波器
F(u,v) =?*(u,v) H(au,av) (12)
进行滤波,上式中
H(au,av) = 4?2a2(u2+v2)exp[-2?2a2 (u2+v2) (13)
vuvuvu,
vuvuvu,
vu
vu
dde)a,a(H)(a}{W
dde)a,a(H)(a}{W
)yx(2i*
y,x,a
)yx(2i*
y,x,a
经滤波后得到 |?(u,v) H*(au,av) |2 (14)
及? *(u,v)?(u,v) H(au,av) H*(au,av) (15)
再经过傅里叶逆变换,分别得到 W{?} 的自相关
a-∞∞|?*(u,v) H(au,av)|2 exp[i2?(ux+vy)]dudv
= Wa,x,y{?}? Wa,x,y{?} (16)
及 W{?}和 W{?} 的互相关
a-∞∞[?*(u,v)H(au,av)][?(u,v)H*(au,av)]
exp[i2?(ux+vy)]dudv
= Wa,x,y{?}? Wa,x,y{?} (17)
其中 W{?} 的自相关给出亮斑 (自相关峰 ),成为目标图像的特征,
由于 小波变换匹配滤波方法 给出经过 边缘增强处理的二维小波变换函数的 自相关,所以有 可能获得更好的识别效果,
由于小波变换匹配滤波器是?*H,所以 H 正是频域中的窗函数,当 a增大 时通带宽度减小,通带的中心频率向低频移动,可以 有效地抑制高频噪声,但鉴别能力下降；反之,当 a 减小 时匹配滤波器将包含更广的频段,
中心频率向高频移动,将 具有更高的鉴别能力,但容易受 噪声的干扰,因为一般的白噪声是与带宽成正比的,
a = l a = 2
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
77
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
78
6.7 光学 Haar小波变换和图形边缘探测
6.7.1 Haar变换
1910年,Haar就首次提出了 Haar变换函数：
h(x)=rect[2(x–1/4)]-rect[2(x–3/4)] (1)
它较易用光学方法实现,因此 Haar小波变换是常用的小波变换之一,h(x)的傅里叶变换为
H(v) = i2?e-i?v [1-cos(?v)]/?v (2)
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
79
6.7.2 Hanr小波变换与边缘探测
Haar小波变换 是信号函数?(x)与 Haar
函数经伸缩后的 母函数 h(x/a)相关的结果,对于一个给定的伸缩因子 a,Haar小波变换的作用如下,
(1) 在小波基元函数 ha,b(x)的正,负半周内对信号进行不加权的积分,这事实上是一个平滑或平均的过程,
(2) 将正,负半周的积分值相减,以上两个作用的综合结果是在平均的意义下求差分,
或求导数,恰恰是测出了图形的边缘,
图 6.17(a)是对一个带有低频噪声的方波 进行
Haar小波变换的结果,
在方波的两个边缘呈现一对峰,极值恰恰指示了边缘的位置,
峰峰实线 为带有噪声的 方波信号 ； 虚线 为 Haar小波变换,
图 6.17(b)是对一个高频噪声的方波 进行
Haar小波变换的结果,
在方波的两个边缘呈现一对峰,极值恰恰指示了边缘的位置,
实线 为带有噪声的 方波信号 ； 虚线 为 Haar小波变换,
峰峰图 6.17(c)方波同时具有低频和高频噪声干扰,
只要伸缩因子 a 选择得当，小波变换仍然有很高的信噪比，峰很尖锐，正确的指示了边缘的所在，充分说明 Haar小波变换的抗干扰能力．
实线 为带有噪声的 方波信号 ； 虚线 为 Haar小波变换,
峰峰第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
83
图 6.18(a)是一个带有噪声的 step函数,g(x)及其 Haar小波变换,噪声是由空间频率为
v1 = 27 及 v2 = 40 的正弦干扰信号及随机本底噪声构成．尽管信噪比很底，小波变换仍然正确地指示了 step函数的边界．
图 6.18 step函数及其 H a a r
小波变换第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
84
图 6,18(b)是信号函数 g(x)和小波函数 h(x)在频域中的行为,实线表示 G(v),虚线表示 H(v)．
G(v)的 主峰主要由 step函数的直流和低频 成分贡献,而 旁瓣 则主要是 step函数 在 边界 的跃变部分贡献的,
零频、低频主峰旁瓣
H(v)频率窗恰恰开在旁瓣干扰信号的主频 v1,v2远在频率窗外第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
85
值得一提的是,两种干扰信号的主频
v1,v2远在频率窗外,对变换没有贡献,而随机本底噪声的谱是宽带的,对结果的影响也不大,
说到底,在空域中干扰信号和随机本底都是,全局,信号,小波变换的局部处理手续恰恰抑制了它们的作用,而突出了图形的局部变化 ——边界的贡献,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
86
6.7.3 二维 Haar小波变换和图形拐角测量二维 Haar变换的母函数定义 如下：
h(x,y) =rect(x - 0.5,y- 05) + rect(x + 0.5,y + 05)
+ rect(x + 0.5,y- 05)+ rect(x- 0.5,y+ 05)
二维 Haar变换母函数如图所示．用 h(x,y)
构作的二维 Haar小波变换，
在探测与 x 轴或与 y 轴平行的边缘时都为零，但测量任意一段既不和 x 轴平行又不和 y 轴平行的边缘时不为零，因而 该变换特别适用于探型图形边缘的
“拐角”,
-1 +1
+1-1
x
y
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
87
图 6.20(a) 给出一个字母,T”作为输入图形,
图 6/20(b) 是由二维 Haar变换测量该输入图形的输出,可以看到拐角测量的效果,
专门进行 x 方向或 y 方向边缘测量的二维 Haar变换,
见图 6.21．
用于测量拐角的 Haar函数常称为,角母函数,,
用于测量边缘的 Haar函数则称为,边母函数,,
根据图形的尺寸,边缘过渡区的宽度及噪声的频谱分布,
可以采用不同的伸缩因子,
以获得最佳的信息提取效果,
图 6.21二维图形边,角测量的 Haar小波变换母函数示意图,(a)角母函数,(b) y边缘母函数,(c) x 边缘母函数,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
89
从一维和二维的边缘探测的讨论中,
我们已看到小波变换的特点,相对于图形整体 而言,边缘 ( 边和 角 )显然 只是局部,Haar小波变换对于图形内外强度不变的或缓慢变化的部分不敏感,即使对一些急剧变化的高频噪声,只要适当选取伸缩因子,Haar小波变换仍有很强的抗干扰能力,但在图形的边缘出现强度跃处,小波变换却出现很强的峰,正确指示边界的位置,可见,在边缘探测这一类问题中,小波变换比傅里叶变换更加适用,
6.7.4 用投影 ——掩模法实现 Haar小波变换本节中介绍一种由 Yang等提出来的在空域中实现
Haar小波变换的方法,投影 ——掩模法,
图 622 用投影 ——掩模法进行 2 维 Haar小波变换
o 为点光源列阵；?l 为小波函数掩模板；
2 为输入信号 SLM；?3 为探测器从点光源 P(x,y)辐射的发散球面波照亮?1上的小波函数掩模板,在 d1的距离上经放大后通过 Ll 投射到
2上,在 Ll 和 L2 间为平行光,没有放大作用,放大率即伸缩因子,a = f /( f – d1 )
因为 P(x,y)一般不在光轴上，使母函数发生平移,
平移量的计算可参见图 6.23．
光线 PQRO2在平面 上的位移量为零，它在 st 平面即?1 平面上的位移量,
bo = O1M + MQ = xd1/f +?/a = x(d1/f + d2/fa ) (5)
式中 a 为伸缩因子,MQ= NR/a =?/a,NR = xd2/f.
bo 放大 a 倍后换成平面即?2 平面上的位移量
b = boa = x(ad1/f + d2/f) = x[d1/(f -d1)+d2/f] (6)
整个光束对于?3 上一点 P’(x’,y’)的贡献为积分
(7)
由 CCD转换成电信号送入数据处理系统,在其中补上因子 a –1/2 即成为小波变换式
(8)
倘若光源是 点光源列阵,我们得到 抽样的小波变换函数 ；若光源是 均匀的连续面光源,我们得到 连续分布的小波变换函数,将掩 模板沿轴向平移,可以改变伸缩因子 a的值,
dd),(a b,a bh)b,b(' yxyx
dd),(a b,a bh
a
1 }{W yx*
b,b,a yx
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
94
6.7.5 偏振编码掩模 ——投影方案中,非相干光 学信息处理用光强大小来表示信号,只能给出 非负的分布,
这显然 违背了小波变换的相容性条件,
若使用 相干光源,可用偏振编码来实现 Haar
变换函数的双极性,
设光源辐射的线偏振光沿 45o方向振动,可以用 Jones矩阵表为
(9)
在掩模中,+1的部分用水平方向的偏振片构成,
用 Jones矩阵表为
(10)
11J? o
00 01P
-1部分用垂直方向的偏振片构成,表为
(11)
输入的光波经过 +1和 -1部分后,分别变成
10 00P
11
11
2
1
P
45
SL MyxJ
J
J
1
0
1
1
10
00
J
P
J
0
1
1
1
00
01
J
P
/4-
o
o
o
表为方向的偏振分析器，可放置沿一前入信号方向的线偏振光。在输和”分别是沿和’
”及
’
掩模板的制作。
编码技术实现了，这样一来，我们就用的相位差为方向的偏振光，它们和分别是沿和
”
及
’
变成振分析器，部分的光波，再经过偏及经过掩模中
H aar
13545J
J
J
1
1
2
1
1
0
11
11
2
1
J
P
J
1
1
2
1
0
1
11
11
2
1
J
P
1-1
oo
-
4/
4/
Wavelet transform and its use
in edge detection
Feijun Song
1 China Daheng Corporation,P.O,
Box 9671,Beijing 100086,China
Suganda Jutamulia
2 Blue Sky Research,3047 Orchard
Parkway,San Jose,CA 95134,USA
ABSTRACT
Application of the wavelet transform to edge
detection is discussed in general,To propose two new
wavelet transforms for edge detection,Modified
Haar’s wavelet transform (MHWT) and Haar-
Gaussian wavelet transform (HGWT).To analyze the
performance of MHWT and HGWT for edge
detection in both space and frequency domains.
Application for lead inspection of surface mount
devices (SMD) in electronic industry by using MHWT
and HGWT is presented in detail.
Keywords,Wavelet transform,Machine vision &
automatic inspection,Edge detection,Lead inspection,
Surface mount devices (SMD)
1,IN TRO DUC TIO N
An i m age g en er ally cont ains s om e slow l y var yin g are as and
sh ar ply var yin g r egi on s,S har p vari eties,si n gula rities o fte n carry
inter esti ng and imp o rta nt i nfo rm atio n.
EDG E D ET ECT ION (I)
Edg es of a n i ma ge a re gene rally local ized by si ngul arities an d i rreg ula r it ies.
Edg e det ecti on is one of the basi c oper ations of al l art ifi cial,com puter,robotic,
or m ach ine visio n syst em s.
Gradi ent est im ations are cl assi cal mo dels which giv e satisf act ory result s in the
case s of sli ghtl y noi sy im ages,
For Noi sy im ages a low -pass fil te r shoul d be inc lud ed,Howe ver,its perf orm anc e
is grea tl y depende nt on the smo othi ng scale fact or of the im pul se respo nse,If the
edge to be d etec te d is shar p,one can choo se a smal l im pul se re sponse,
EDG E D ET ECT IO N (II)
Bu t if th e ed ge is sm oo th,a la rge im pu lse re spo nse will gi ve m ore
sat isfa ct ory re sul ts.
Th er ef ore,sin gl e sm oo th in g scale of im pu lse re spo nse is oft en
in suff ic ie nt in m an y ap pl ica tio ns,an d m ul ti-sc al e or va ria ble im pu lse
resp on se fo r e dg e d et ec ti on is ve ry at trac ti ve.
WA VEL ET TR A NSF O RM
M all at [ 1] su g ge ste d t hat th e w a vel et tr an sf orm is very c l os ely rel a te d
to a m ult i- s cal e i mp uls e r es po ns e for e d ge d etec ti o n,
An at tr acti ve f eatu r e of t he w a vel et tra ns fo rm is th at it ca n wel l l oc ali z e
th e re gu la rit y of sign als in b oth spa c e an d f req u en cy d om ain s,By
va ryin g it s dil ati on p ar a met e r,th e wa v el et tr an sf orm a cc o rdin gl y
ch an g es i ts s mo oth in g sc ale f act or au to m ati c all y.
2,NE W W AV E LE T T RA NS FO RM S
H aa r’ s wa ve l et is co ns t ru cte d by t wo re ct an gu l ar fu nc t i on s as ex pr es se d i n
Eq,(1) an d sh ow n in Fi g,1,
4
3
2
4
1
2 xr e c txr e c txh
,(1 )
wh er e r ect ( x ) i s defi n ed a s
2/1||0
2/1||1
x
x
x
H AA R W A VE L E T
Its Fou rier tran s for m i s
)2(,
c o s1
2
i
ieH
,
It h as a ze ro dc c om po n ent,
Fi g.1,H aa r’ s wa ve let (a ) a nd it s s pe ct ru m (b )
x
h ( x )
1
( b)
H ( )?
.
( a)
TWO NEW WAVELE TS
We pr opos ed tw o new wav ete s,Mod ified Haar ’ s wavel et and Haar-Gaussian
wav el et,[4,5],d efine d as:
)0(' qs
s
qx
rec t
s
qx
rec txh
s
,(3)
)0(e x pe x p"
22
qs
s
qx
s
qx
xh
s
(4)
q,off set param e t er s,di l at i on fact or
Modifi ed H aar’s w av elet an d Haar -Gaus si an wavel et are b ot h anti-symme tri c as
show n i n Fig,2.
Wit h an ad ditional of fset
para met er q,m odifie d Haa r’ s
wav el et and H aar-Gau ssian
wav el et ar e ver y powe rful in
noisy edg e dete ction and
featur e ext raction,
Fig,2,Modifi e d Haar ’s wave le t
( sol id line ) and
Haar -Gaus sia n wave let
( das hed line)
- q q
3,NOISY EDGE DE TECTIO N
Math em at ic al mo de l of edge functi on is expr essed i n Eq,(5)
2
1t a n
1
x
xf
,(5)
Wi dth of ed ge? is
given by
,( 6)
The ed ge functi on is
depict ed in Fi g.3.
Fig,3,Model of an edg e
with w id th?,
NOISY EDGE FUN CTION (I)
No is y edge functi on is furt her s i m ulated as f ollow s.
r n dxx
x
xf
2211
1
2s i n2s i n
2
1t a n
,(7)
Last t hree ter ms repr esent the sinu soid int erferen ce and rando m n ois e,
Si mu lat ed nois y edge functi on is depict ed i n Fi gs,4(a),(b),(c) an d ( d),t ogether with
the H aar-Gaussi an wavelets of various d il ations s.
The vertical li nes ( da shed) specif y t he sp at i al win do w s of the wavelet,
( a) ( b)
(a) s=0.63 (b) s=0.42
NOISY E DGE FUNCTION (II)
Fig,4,No is y e dge functi on a nd Haar-Gaus sian w avelets with various dil at ions,
(a) s = 0,63,2s = 9.5? ; ( b) s = 0,42,2s = 6,4? ;
(c) s = 0,10,2s = 1.6? ; (d ) s = 0.053,2s = 0,8?
2s,w idt h of th e spa ti al window, is de fi ned as follow s.
2/1
2
","
","
2
xhxh
xhxxh
s?
,(8)
where ( f,g ) repres en t s the inn er product o f f un cti on s f and g,
(c) s=0.01
(d) s=0.05 3
HAAR-GAU SSIAN WA VELE T TRA NSFORM OF TH E
NOISY EDGE FUN CTION (I)
Ha ar- G aussi an wa vel et tr ansf or m ( HGW T ) of t he noi s y edg e
functi on f ( x ) is a cor r elati on of t he sign al wit h the wa vel et:
d
s
x
hfxhxfxfW
ss?
)("
*
,
(9)
HGW T for diff er ent di lat i ons s (wh i c h is equi val ent t o th e
smo oth in g sc ale of im pul s e r es pons e dis c uss ed pr e vio usl y) ar e
depi ct e d in F ig,5.
M axi ma of th e W T in di cat e
th e locati on of the e dge,
Fi g,5,HG WT of the no i sy edge
fo r di ff erent dil a ti on s,
(a) s = 0.63 (cu rv e wit h t he
br oad est peak),
(b) s = 0.42 (cu rv e wit h t he l ess
br oad peak),
(c) s = 0.10 (cu rv e wit h nar ro w
peak and sm al l ri pp l e s),and
(d) s = 0.05 3 (cu rv e wit h sev ere
ri pp l es),
HAA R-GA USSIAN WAVELET TRANS FORM OF THE
NOI SY ED GE FUNCTI ON (II)
Charac t eris tic d i me nsi o n o f the noise is gi ven as:
2.0/1
,(10)
,th e freq ue nc y o f the sinu soi d al n o ise.
,
s,d ilat io n of WT
2s,wid th of sp atial wi n d o ws
,char acteriz ed n o is e di mens io n,
When th e wid t h o f the spa tial w i n do w is clo se t o th e ch aracteristic
d imension of t he n o ise multiple ma xim a a re found in the wavelet
transfor ms and the fal se response w o ul d oc cur.
4,FREQUEN CY DOMAIN ANALYSI S
The edge detection based on the wa velet transf orm in space domain can
be analyzed in frequency domain as we ll,The frequency spectra of both
the noisy edge function and the Ha ar-Gaussian wa velet of differ ent
dilations s are depicted in F ig.6 (a) to (d).
S 2s
0.63 3.0
0.42 2.0
0.10 0.50
0.053 0.25
0.2
Tab le 1
Fig.6,Spectra of noisy edge functi on (sol id line),H aar-G aussian
wavelet (dashed li ne) and frequency fi lt ering window (Verti cal li nes),
(c) s = 0.12
(a) s = 0.63
(d) s = 0.053
(b) s =0.42
FREQUENCY DOM A IN ANAL YSIS ( I)
Ha ar- Gaussi a n
wa vel et acts a s a
low fr e quen c y
ban d- pass fi lt er
wi th no dc
comp onent i n
fre que ncy dom ain,
The s pect r um of the r a mp edg e fu nct io n consist s of a hi gh dc
comp on ent a nd a se ri e s of low f re quen cy l obe s,The s harp est vari ati on
of t he il lu mi nant cur ve (t he e dg e) c ontri but es mo stl y t o th e l ow
fre que ncy l obes,
The ba nd-pa ss fi lt er of th e Haa r -Gau s s ia n wav e let block s th e dc
comp on ent but pas s es th e low fr equ enc y l obes belo ng t o th e e dge
sign al a nd th us e nh an ces and ext ract s th e e dge tr ansi ti on
in form ati on.
FR EQ UEN CY DO M A IN ANAL Y SIS (II )
Q q
FR EQ UEN CY DO M A IN ANAL Y SIS (II I)
In Fi gs,6( a) a nd (b),th e ban d-pas s fil t er of th e HGW T blocks th e tw o hi gh-
fre quen cy pea ks of th e sin usoi d noi se fun cti ons,These tw o pe aks a re fa r
bey ond t he cut -of f f r equ ency of t he ban d-pass fil ter,
Al so,t he spe ctru m of ra ndom noi se is broa dl y spr eadin g,havi ng le ss in flu enc e
to t he wave let tr a nsf or m.
(a) s = 0.63
(b) s=0,42
In co nt ra st t o Fi g.6 (a ) an d (b ),t he fr eq ue nc y
wi nd ow of Fi g.6 (d ) co nt ai ns t he no i se
ch ar ac t er i st i c fr eq ue nc y,whi ch pr od uc es t he
m ul ti- m ax im a i n t he sp ac e do m ai n an d m i gh t
re su l t i n som e fa ul t l oc ati ons i n ed ge
de t ec t i on,
Fi g.6 (c ) re pr es en t s a cr i tic al
si t uati on whe re t he pe ak s j us t l oc at e at
t he m ar gin of t he fr eq ue nc y wi nd ow,
Th ei r co nt ri but i on i s i nsuf fi c i en t t o
i m pa ct th e ed ge det ec tio n.
Fi g,6 (c)
Fi g 6 (c )
Fi g,6 (d)
FR EQ UN CY DO M AIN AN AL Y SIS ( IV)
CO NC L UT IO N
The fr eque n c y dom ai n a nal ysis i s i n a g o od a gre e m ent wi th t h e spa c e d om ain a n alysi s.
In t he fina l a nal ysis th e r an d o m an d sin u soid al n oise h av e in fin it e e xt en sion s in sp ace,
but th e wa v el et tr an sfo rm loc ali z es th e si gn al in b oth s pa c e an d fr equ en cy dom a i n s,It
is th ere for e i ns e nsit iv e t o the n ois e a nd abl e t o accur a tel y in di cat e th e lo c ati on of th e
edg e.
5,AP PL I CA T IO N T O LE AD INS P EC T IO N O F
SU RF AC E M O UN T D EV I CE S
BA CK G RO UND
In t he last d e ca de th e v olu me p ro d ucti o n of s urfa c e mo u nt de vic e ( S M D ) has gr o w n
trem e nd ou sly a n d t h er e is a t re n d t o r epl ac e t h e th ro u gh -hole c o m po n e nts b y S M D.
T HR OU GH- HOL E C OM P ON E N T S A N D
SUR F A C E M OU N T D E V I C E ( SM D )
S u r fac e m ou n t d e v i c e (S M D ) T h r o u gh -h o l e d e v i c e
Th i s ne w t e c hn ol og y p r od uc e s b e ne f i t s i n
M i n i at u r i z at i o n
M u c h l ar ge r i n t e gr at i o n
c os t r e d u c t i o n,a nd
d e s i g n f l e xi b i l i t y,
COP L ANA R I T Y,P I T CH AND F OO T ANGL E OF S M D
H ow e ve r,t he f ol l ow i ng f a c t o r s of S M D h a ve put m uc h s t r e ss on t he p r oduc t i on
t e c hnol ogy,v e ry f i n e r pi t c h ( <0.5 m m ) a nd m u c h h i g h e r pac ka ge l e ad c o u n t s
( m o r e t h a n s o m e hund r e ds l e a d s pe r c h i p ),
To a s s u r e r e l i a bl e s ol de r c ont a c t s,i t i s r e q u i r e d t h a t t he l e a d c opl a n a ri t y a nd
pi t c h s houl d b e 1 00 m i c r ons,The f oot a n g l e s houl d b e l e s s t h a n 10 de g.
The I C m a nuf a c t u r e r s a r e r e s pons i bl e f o r a c c u r a t e l y i n s pe c t i ng t he 3- D
pos i t i o n s of t h e l e a ds f o r a l l ou t goi ng c o m p one nt s,
C
P
C o p l an ar i ty ( C ),p i t c h ( P ) an d fo o t an gl e (? ) o f S M D
DUAL SH A DO W TE CH NIQ UE F O R LE AD INSP ECT IO N
S everal te chni qu es h ave b een pr o pos ed f or th e le ad i ns pect io n o n th e pr o du ctio n
line,T h e du al sh ad o w te chni qu e [6 ] is cos t-effec t i ve,r elia ble,an d acc ept ed b y
most m anuf act ur ers,
A s che matic di agr am of th e lead i n s pect ion syst em is sho w n in Fi g,7.
An S M D is p ut on th e ped estal,All its l ead s do n ot t ou ch th e to p pl an e of a
defu sin g scr een s o as t o m ake a n on des tru ctiv e t estin g sy ste m,
The col lim ate d beam fir st illu min ate s the l ead s o f on e sid e from l ef t,and then
the se cond col limat ed beam illumi n ates th e sa me lead s from r i gh t,
Tw o sh ad o ws ar e pr o je cte d on t he su rfa ce of a de fu sin g scr een.
The s hado ws ar e th en im ag ed by a le ns o nto a CCD an d grabbe d and
pr oces s ed b y a d ata pro ces sin g sys tem,
L E A D I N SP E C T I ON SYST E M
F i g,7,S c h e m at i c d i ag r am o f th e l e ad i n s p e c t i o n o f S M D
b as e d o n d u al s h ad o w te c h n i q u e,
c ol l i m a t e d be a m s c ol l i m a t e d be a m s
l e f t pr o j e c t e d i m a ge
d e f u s i ng
ed
s c r e e n
CCD
l e ns
p e de s t a l
s c a nn i n g
l i ne of C C D
r i g h t pr o j e c t e d i m a g e
B ASI C M E ASUR E M E NT OF L E AD I NSP E CT I ON
A b as i c m e as u re m e n t of t h e i n s pe c t i on is t o ac c u r at e l y de t e rm i n e t h e
pos i t i o n s of t h e c e n t e rs of l e ads i n a s c a nn i ng l i ne of C C D,
A n e d ge de t e c t i o n t e c hni q ue i s g e ne r a l l y a ppl i e d t o de t e r m i ne t he f ro n t a l
an d t ra i l i n g e dge s of t he s h a dow of a l e a d,The c e n t e r po si t i o n c a n t he n b e
de t e r m i ne d b y t he pos i t i on s of e dge s.
NUMERI CAL DIF ERENT IAT ION A ND W AVELE T TRA NSFORM
Numerical dif ferentiation is often unab le to use for the following reasons:
1,The edge is a ramp that canno t be d efi ned definitel y,
2,T he detect ed gray scale c urve o f leads is corrupted by nois e,
3,T he cu r ve flu ctuates du e t o t he aging of the li gh t source an d especiall y du e t o t he
va ri ation of the backgro und brig htness.
To overcom e t he pro blem,w e ap ply the newly prop osed MHWT and H GWT,T he
image of a lead,wh i ch likes a va l l ey in the gray s cale cur ve (ima ge s of le a ds )
shown i n Fig,8,is al m ost s ymmetric,
Therefore w e do not de t ec t the edges of t he valleys.
Ins t ead,w e de t ermine t he symm etri ca l center s of the valleys dire ct l y using the
MHWT and HG WT,The algorithm is called sym me t ri c evaluati on,
s
qq
Fig,(8)
O FF SE T A ND DIL AT IO N P AR AM ET ERS
In t he WT,th e o ff set p ar ame te r q re p r es e nts t he offs e t of th e e d ge t o t h e
sym m etri c al cent e r,a nd w e s ele ct t he dil a ti on pa ra m eter s to ma t c h th e wi dt h of
the ra m p o f edg e ap pr oxi mat ely,
q q 2s
AD AP TI VE L YL Y S M OO TH I NG S CA LI NG
Th e wa ve let tra ns fo rm pr od uc es m ul ti-s cale im pu lse res po ns e.
W e ha ve fu rth erm or e de sig ne d a sp ecia l m eri t fu nc tio n to
ev alu ate th e in sp ecti on pr oc ess in g,W he n th e wi dt h of th e
ed ge is va rie d du e to so m e re aso ns,th e di la tio n pa ra m eter s is
ch an ge d to m at ch th e wid th of th e ed ge ad ap tive ly,ca lled
ad ap tive ly sm oo th in g s ca lin g.
EXP ERIM ENT RES ULT
Fi g,8,The de te ct ed gr ay sc al e cu rv e of th e l ea d sha do w (l i ne
wi th + si gn s) an d i ts m od i fie d H aa r wav el et tr an sf or m
(dashed line) and H aa r- G au ssi an wav el et tr an sf or m (s oli d
l i ne ).
Bec a us e b oth mo difi ed Ha ar ’s an d Haar - Gau ss ian
wav el ets a r e an ti -sy m m etri c,a n d t h e vall e y of th e si gn al is
alm ost sym m etri c t h e wa v elet tr a nsf or ms of t h e gra y sc ale
curv e ( c orr elati on o f si gn al a n d th e w a vel et ) pro d u ce
dire ctl y th e z er o-c ro ssi n gs at t he c ente rs o f th e v a ll e ys
(i ma g e of l eads ),
T he m a xim a an d mi n im a o f th e w avel et t ran sfor ms
re pre sen t th e e d ges of th e lea d sh ad ows,an d t he zero-
cr ossi n g s accu r atel y in di cate th e c en ter s of l eads,
A P P L I C A T I O N O F W A V E L E T T R A N S F O R M T O
L E A D I N S P E C T O N O F S M D
T he w a ve l e t t r a ns f or m s a r e r e p r e s e nt e d b y
T he w a ve l e t t r a ns f or m i s a c o n t i n u o u s f u n c t i o n o f t h e s p a c e v a r i a b l e
x,a nd i s a d i s c r e t e f u n c t i o n o f t h e f r e q u e n c y v a r i a b l e ( d i l a t i o n p a r a m e t e r ) s,
A s p e c i a l m e r i t f u n c t i o n i s d e s i n e d t o o p t i m i z e t h e d i l a t i o n p a r a m e t e r s
a s di s c u s s e d b e f or e,
d
s
x
hfxhxfxfW
ss?
)("
*
n
ssss,...,,
21
P H O T O G R A P H Y O F C O N S T R U C T E D S Y S T E M F O R L E A D
IN S P E C T IO N
T h e p h o to g ra p h o f t h e r e al se tu p i s s h o w n i n F ig,9,
T h e i n st ru m e nts h a v e b e e n p ut in to m a n u fa ctu ri n g,
CO NCL U DING RE M A RK S
W e ha ve pro po sed ne w m od ifie d Ha ar ’s wa ve le t tra nsfo rm an d
Ha ar -Ga ussi an wa ve le t t ra nsfo rm,
Th e pro po sed wa ve le t tra nsfo rm s ca n be ap pli ed to th e ed ge
de te cti on,
Bo th wa ve le t tra nsfo rm s ar e al so pro ve d ve ry po we rfu l in find in g
th e sym m etri c ce ntra l po siti on of no isy ba rs,whi ch can be ap pli ed
to th e r eal le ad insp ec tio n o f su rfa ce m ou nt devi ce s (SM D).
Th e in stru m ents ha ve be en put in to m an ufa ctur in g an d ap plie d in
th e p rod uc ti on lin es.
DIF FER ENCES BE TWEEN FOUR IER TRAN SF ORM AN D W A VEL ET TRAN SFOR M
T he Fouri e r transform (FT) e xtrac ts the p e riod ic signals on infinite duration,while t he
wa velet transform (WT ) look s a t the window specified by its dilation para m e t e r s a nd
c a n be shift e d by its transla tion para m e ter x,T h e FT i s suitab le for proce ssing
signals on infinite or long r ange duration,but the WT is c a p a ble of proce ssing
transition o r short -term signals.
For a c tual a ppli cation the a ppare nt non- z e ro c o m ponents of a wa v e let only e xist in a
very sm a ll region in both spac e a nd freque ncy dom a ins,An a ttrac tive feature of the
WT is that it lo c a liz e s the proce ssing of signal in both spac e a nd freque ncy dom a ins.
T here fore WT is a space - and freque ncy localizati on operat or,
DIFFE REN CES BET WEE N FO UR IER TRA NS FO RM AN D W A VELET
TRA NS FOR M
The fr eq u e n c y v aria bl e? a n d s p ace vari a ble x ap pe a r i n th e F T a nd
inv erse F T,re sp ec ti ve l y
.)2e x p ()()(
,)2e x p ()()(
dxjFxf
dxxjxfF
The WT p r o vid es in si gh t th at simu lt ane ousl y com bin e s th e
vari a ble s in b oth sp a c e and fr eq u ency dom ain s,i.e,th e
tran slati on p aram eter x a n d d il atio n p aram eter s ap pe a r i n th e WT
simult an e o u sly,
.)()}({
*
d
s
x
hfxfW
s
第六章光 学 小 波 变 换第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
2
第六章 光学小波变换
6,1 引言
6,2 从短时傅里叶变换到小波变换
6,3 小波变换的定义和性质
6,4 实现一维小波变换的光学系统
6,5 用多通道匹配滤波实现二维小波变换
6,6 光学小波变换匹配滤波器在图像识别中的应用
6,7 光学 Haar 小波变换和图形边缘探测
6.1 引 言如果 g(x)是一个时域或空域中分布在 (-?,?)中的稳恒过程或稳定分布,则傅里叶分析给出了近乎完美的结果,然而,在自然界和科学技术中还有大量信号,它们具有局部的或定域的特性,例如语言信号,声纳信号,
各种电脉冲等 。 这些 信号只出现在一个短暂的时间问隔内,此后很快衰减到零,称快速过程或暂态过程,一个很短暂的信号,可以 称为,小波,信号对于局部信号或暂态过程,傅里叶分析就不完全适用,首先,我们仅对?t内的时间信号感兴趣,没有必要在过去,现在及未来的无限长时间范围内对信号进行分析；类似地,在处理定域于?x 内的空间图像时,
也没有必要对全平面内的信号进行全面的分析,
在许多情况下,在?t 或?x 以外的信号是未知的,
它可能是零,也可能是背景噪声；对它们我们不太了解,测不准,或不感兴趣,如不加选择地把 (-?,?)内全部信号进行傅里叶处理,还可能产生较大的误差甚至错误,此外,一个局部的信号在?t 或?x 以外较远处几乎完全等于零,当用它们的频谱来恢复或重构这些信号时,在?t 或?x 外很远处也会出现一些非零的分量,它们一般不是信号,而是 在傅里叶逆变换中频域综合不够充分而产生的噪声,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
5
在一些课题中,我们往往不满足于了解信号在全部区间内的综合的频谱分布,而 希望了解某一区间或某些区间内信号对应的频谱,例如在地震勘探中,为了分辨分层的地层和矿床结构,我们需要在时域和频域中仔细分析不同时刻的信号在不同频谱区间中的行为,而傅里叶分析只能提供在长时间内的信号整体的频谱,显然不能满足我们的要求,
近年来发展起来的小波分析,正好克服了傅里叶分析的上述缺点,它和傅里叶分析的一个重要区别,在于它恰恰适用于处理局部或暂态信号,因此,小波分析成为信号分析,图像处理,
数据压缩,语音信号分析等领域中的重要工具；
在地震勘探信号处理,边缘探测,语音信号合成中则有特殊的用途,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
6
6.2 从短时傅里叶变换到小波变换
6.2.1 短时傅里叶变换 (STFT)
为了有效地提取一个局部信号 g(x)的信息,引入一个局部化的变换,
所谓局部化,包含两个要素：
第一,被分析的区间要有一定的宽度
x,我们仅对?x及其附近的信息进行处理；
第二,被分析的区间有一个中心坐标 xc，
当 xc 改变时,就可以提取不同的信息,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
7
为了实现局部化,一个有效的方案是在傅里叶变换中加一个 窗函数 w(x)：
Gw(v,xo) =?∞- ∞g(x)exp(-i2?vx) w(x- xo) dx (1)
= [W(v)exp(-i2?vxo)]* G(v) (2)
式中 W 和 G 分别是 w 和 g 的傅里叶变换,
只要 w(x)和 W(v)有足够快的衰减速度,窗函数就是一个局部化的函数,
函数 f (x) 和 g(x) 的 内积的定义,
( f (x),g(x) ) =?∞-∞ f (x) g(x)dx
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
8
窗函数的中心 xc定义,
xc = (w(x),x w(x))/(w(x),w(x))
=?∞- ∞ w(x) xw(x)dx /?∞- ∞ w(x)w(x)dx
式中 (*,*)表示函数 f 和 g 的内积,xc与 xo不一定相等。
窗函数的宽度 则定义：
w = 2[(w(x),(x- xc )2 w(x))/(w(x),w(x))]1/2
= 2[?∞-∞w(x)(x-xc )2w(x)dx /?∞-∞ w(x) w(x)dx]1/2
注意,它是节 1.2.4中所定义的信号空域宽度的两倍 。
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
9
Gw(v,xo) =?∞- ∞g(x)exp(-i2?vx) w(x- xo) dx (1)
由于窗函数具有局部处理的功能,因此 (1)式定义的变换称为 短时傅里叶变换 (STFT)．
STFT和 FT的一个重要区别,
频率变量 v 和坐标变量 xo同时出现在变换函数中,
在 STFT中,窗口宽度则隐含于 Gw(v,xo) 内,
正是 xo和窗口宽度?w,使这一变换具有局部处理的功能,改变 xo,窗口就在空域中移动,以获取不同区域的信息,xo 通常称为 位移因子 ；?w 则限制了被处理空间的范围,
频率窗中心,vc = (W(v),vW(v))/(W(v),W(v)) (5)
频率窗宽度,
W=2[(W(v),(v- vc)2W(v)w(x))/(W(v),W(v))]1/2 (6)
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
10
当?w和?W 都有限时，我们称函数 w(x)在空域和频域同时局部化．
w?W称为 空间 -频率窗,它限制了空域和频域中被处理区域的范围。根据?w,?W的定义及测不准关系式，并注意信号宽度定义的区别，我们就有,?w?W≥1/? (7)
当 高斯函数取为窗函数 时,(7)式中的等式成立,这种情况下 STFT具有最小处理区域,
STFT的局部性,其特征在于处理过程限制在空间 -频率窗内进行,且窗的位置是可变的,然而无论?w还是?W都是常数,不会随信号中心频率的变化而变化,这使 STFT在处理一些奇异性的信号时显得无力,恰恰是在这一点上,小波变换 具备比 STFT更强的功能,
6.2.2 Gabor 变换早在 1946年,Gabor就提出了下面的变换
(8)
称 Gabor变换,其中?和 b为变换的参数,上式又可表为
(9)
式中
(10)
因此 Gabor变换就是高斯窗短时傅里叶变换,
窗函数中心坐标 xc = 0,窗的宽度?w = 1.414?
w?(x)的傅里叶变换,W?(v) = exp(-2?2?2v2) (13)
也是高斯函数,频率窗宽度?W = 1/ 1.414 (14)
因此有?w?W = 1 /?,
图中将空域和频域同时表达出来,称 空间 -频率坐标系,空 -频窗 则表示为图中的一个矩形,Gabor变换空 -频窗的高度和宽度都是恒定的 。
Gabor变换在频域中的表达式,[式中? = 1/(2)2 ]
可见 Gabor变换在频域和空域中的表达式具有相似的形式．
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
13
Gabor变换的特点,
(1)实现空域和频域处理的局部化 [中心位于 (b,v)
空 -频窗为 1.414?× (1/ 1.414) = 1/? ]
(2) 由 (8)式和 (15)式,
可看出变换是参数?,b 和变量 v 的函数。上式给出的 积分是一个调制包络,载波 exp(-i2?vb)的频率 (即中心频率 ) v 与参数?无关,不会随? 的变化而变化，这正是所有短时傅里叶变换共同的 缺点 。
(17)
6.2.3 Morlet小波变换为了克服 Gabor变换中窗口尺寸不能变动的缺点,
Gabor变换的 基元函数 可改写为子函数定义信号函数 g(x)的 Morlet小波变换,
称 变换的母函数,
引入参数 a,b，
Morlet小波变换与 Gabor变换的实质性差别,
小波变换,vc= v/a,?w =1.414?a,?W=1/1.414a,
当 vc 增高时 (a 减小 ),?w 变小而?W增大，可处理更多的高频信息；当 vc 降低时 (a 增大 ),?W变小而
w加宽，可容纳足够多个空间周期，以保证处理精度．
Morlet小波变换,在处理低频信号时空间窗自动加宽,
在空间窗范围内包含的信号空间周期相同,这就保证了小波变换以同样的精度去处理不同中心频率的信号,这正是小波变换与短时傅里叶变换的根本区别,
Gabor变换,窗的宽度是常数，当 vc增高时，一定宽度的空间窗内包含的空间周期增加，所以 变换的精度是随频率而变化的 ；
小波变换的空间-
频率窗第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
16
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
17
6.3 小波变换的定义和性质
6.3.1 小波变换的定义母函数 h(x)的基本小波函数 ha,b(x)定义为
(1)
式中 b称为小波变换的 位移因子,a > 0 称为 伸缩因子,上式表明基本小波是母函数经平移和缩放的结果,基本小波又简称 小波,
信号函数 g(x)的 小波变换定义 为,(2)
由于相关运算较易用光学相关器进行,因此小波变换可以用我们已熟悉的光学相关系统来实现,
)a bx(ha1)x(h b,a
)b(g)ab(ha1dx)x(g)a bx(ha1)x(g),x(h)}x(g{W *b,ab,a
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
18
Morlet小波的 母函数 是子函数 是
Mor1et小波的基函数式中,m,n = 0,1,2,… ； a = aom ； b = nbo.
)2 te x p ()t2ie x p (2 1)t(h 2
2
0
)a bt(ha1)t(h n,m
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
19
当 m = n = 0时,a = 1,b = 0,即为母函数；
当 m = 1,n = 0时,a = a0,b = 0,对应的一阶小波函数为当 m = -1,n = 0时,a = a0-1,b = 0,有负一阶小波函数
h(t),h1,0(t),和 h-1,0(t)已分别在上图中画出,
)a2 te x p ()a t2ie x p (
2
1
a
1)t(h
2
0
2
2
0
0
0
0,1
)2 tae x p ()ta2ie x p (
2
1
a/1
1)t(h
2
22
0
00
0
0,1
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
20
并非任何函数都可以作为小波变换的函数 h(x),h(x)必须在?x时衰减到零,实际使用的小波变换母函数 h(x),当
x 时迅速衰减,使它的不显著为零的分量只存在于一个很小的区间内,这正是,小波,名称的来由 。 实际上,也只有迅速衰减的小波才使变换 (2)式具备局部化的特征,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
21
6.3.2 小波变换在频域中的表达式在频域中,小波可表示为式中 H(v)是小波母函数 h(x)的傅里叶变换,在空域的扩大 x/a 等价于频域的压缩 av,空域的位移 b 等价于频域的位相移动 exp(-i2?vb)．
由 Parseval定理,小波变换的傅里叶变换为
(4)
其中 H 和 G 分别是 h 和 g 的傅里叶变换,上式表明 信号 g(x)的小波变换可以用 4f系统实现,
)a(H)b2ie x p (adx)x(h)x2ie x p ()x(H b,ab,a
d)(G)b2ie x p ()a(Ha)}x(g{W *b,a?
6.3.3 逆变换和相容性条件小波变换 (2)式的逆变换定义为其中 Ch满足条件
(6)式称,相容性条件,,是逆变换存在的条件,
下面我们来证明它,ha,b(x) 可表为以 (4),(7)式代入 (5)式,得到上式中最内部的积分为?(v’ - v),因此上式成立的条件是即 (6)式．当 v = 0时,相容性条件 要求,H(0) = 0,(10)
即小波函数没有零频分量．由于 H(0)=?-h(x)dx = 0 (11)
意味着 h(x)必须是振荡函数,平均值为零,其傅里叶谱直流分量为零,
6.3.4 正则性从理论上讲,任何满足相容性条件的函数都可当作小波变换的母函数,然而在实用中,为了使变换具备局部化的功能,h(x)和 H(v)在空域和额域中都是迅速衰减的,
它们不显著为零的分量分别分布于空域和领域中的原点附近,
此外，要求 Wa,b作为 a的函数，应当是充分光滑的，
当 a?0 时,Wa,b?0，即要求 Wa,b在 a=0附近是正则的,设
b = 0，则有将 g(x)在 x =0的邻域内展开成泰勒级数：
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
25
代入 (12)式得到式中此是小波函数 h 的 n 阶矩：
Mn =?- h(?)?nd? (n = 0,1,…) (15)
由 (11)式 Mo = H(0) = 0 (16)
设 Mp = 0 ( p = 0,1,…,n) (17)
则在 0 的邻域内第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
26
随着 a? 0,Wa.o? 0 的速率为即对于一个足够平滑的函数 g(x),Wa.o{g}以 a n+1/2
的速率随 a 趋近于零，称它为 n 阶小波函数,
由节 1.1(24)式,Mp = 0 ( p = 0,1,…,n)意味着 H(p) (v) = 0 ( p = 0,1,…,n)
表明 v = 0是 H(v) 的 n 阶零点,
此外,(16)式 (即相容性条件 )保证 Wa.o随 a 趋于零的速度的下限为 a1/2．
6.3.5 小波变换的空间 -频率窗和处理过程的局部化小波变换在空域中的处理局限于空间窗内
(23)
小波变换在频域中的处理局限于频率窗内
(29)
空间 -频率的处理就局限于空间 -频率窗内：
(31)
小波变换处理过程的特点：
(1)空间窗宽度 a?w 和频率窗宽度?W/a 均随 a 的变化而变化,窗的面积?w?W 与 a 无关,
(2) 中心频率 v /a 与带宽 (即频率窗宽 )之比
Q = (v /a)/ (?W/a) = v /?W (32)
与中心频率大小无关,仅取决于 H(v),Q是测量精度的特征量,上式表明 小波变换的测量精度与频率无关,
当 v/a增大时 (a 减小 ),频率窗自动变宽,使小波变换在不同频率下具有相同的检测精度；反之,当 v/a 减小时 (a 增大 ),空间窗自动加宽,以容纳同样数目的信号空间周期,有人把小波变换的这种性能比喻为,自动变焦,(zooming)，
伸缩因子 a 常称小波变换的 频率变量,位移因子 b 则称为 坐标变量,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
29
6.3.6 常用小波函数
1,Haar小波
Haar小波是双极性阶跃函数它以 t = 1/2为中心的奇对称实函数,满足小波存在条件,
)]43t(2[r e c t)]41t(2[r e c t)t(h
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
30
Haar小波的傅里叶变换 是它的模是正的偶函数，以 v = 0 为对称轴,相位因子 exp(-i?v)是由于 h(t)是以 t =1/2为中心的奇对称性所引起的．
Haar小波对于离散的缩放因子和位移量是正交的，其傅里叶谱的振幅 |H(v)|随 1/v 很慢地收敛到零．
/)]c os (1)[ie xp (i2)(H vvvv
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
31
2,Morlet小波
Morlet小波是由分析声像技术引人的,
可表示为其实部是余弦 -高斯函数,
Morlet小波的傅里叶谱 是平移到 vo 及 –vo
处的两个高斯函数，即显然,H(v)是正的实偶函数．
)2te x p ()t2ie x p ()t(h
2
0
]})(2e x p [])(2{ e xp [2)(H 202202
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
32
图,h(t)的实部和 h(t) 的傅里叶谱由于 H(0)≧ 0,所以小波变换的存在条件没有满足,但是对于较大的 vo,H(0)十分接近于 0。 在数值计算时可近似地看作为 0.
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
33
3,Mexican – hat 小波这种型如墨西哥帽的小波函数被广泛用于零交叉多分辨边缘检测,其母函数实际上是高斯函数的二阶导数,即该函数是实偶函数，满足小波变换存在条件,
Mexican - hat小波的傅里叶变换是它也是实偶函数,
)
2
t
e x p ()t1()t(h
2
2
)2e x p (4)(H 222
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
34
Mexican - hat小波函数及其傅里叶变换高斯函数的高阶导数也可以用来作为小波函数,高斯函数 n 阶导数的傅里叶谱 H(v)
将是乘以 ( i 2?v )n 的高斯函数,所以 H(0) =
0，满足子波变换的存在条件。 Mexican –
hat 小波变换收敛速度很快．
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
35
4,Meyer 小波
Meyer小波是用其傅里叶变换来定义的,即
H(v) = exp(-i?v )sin[ω(v)]
式中,ω(v)是偶对称函数，如图所示．
图中,AB弧段有一个对称中心，在 v =1/2 和
ω(1/2) =?/4 处,即
3231),(2)1(
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
36
而 BC弧段具有相同的形状，是 AB弧段的翻转和展开，其对称中心在 v = 1 和 ω(1) = π/4，即如果 H(v)中相移因子 ex(-i?v)被忽略，则 H(v)为实偶函数．相移因子的作用，是使 h(t)在时间轴上的位移为 t =1/2．
Meyer小波可表示 为它也是一个以 t =1/2 为对称点的实函数，衰减较快,Meyer已经证明，具有离散缩放和位移因子的这种小波，可以构成一组正交基底．
3231),(2)2(
d])21t(2c o s [)](s i n [2)t(h
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
37
6.3.7 离散小波变换 (DWT)
正如傅里叶变换一样,小波变换也可以分为两种形式：连续的和离散的变换,
对于 离散小波变换,我们将假定信号 g(x)也是离散的,离散小波变换的形式为对应 离散小波逆变换 为
a/]a/)xx[(h)x(gW
x
jjiij
a/]a/)xx[(hW)x(s
i j
jjiij
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
38
6.4 实现一维小波变换的光学系统从小波变换的定义可知,N维信号函数的小波变换是 2N维函数,因此计算的工作量很大,尽管专门用于小波变换的超大规模集成电路已有报道,但人们仍在考虑用光学系统来实现小波变换,
因为光学信息处理器具有高度的并行处理性能,
一维小波变换光学处理系统第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
39
xy-L2-uv 构成 x 方向的一维傅里叶变换器。在
SLM1上输入信号 g(x)，经 L2 的变换作用，在 uv
平面上形成它的傅里叶谱 G(u)．
在 uv 平面上放置第二个空间光调制器 SLM2，
它被分成 M个沿 u方向的带状区域，这些带状区域中分别显示具有不同伸缩因子 am 的基元函数 h
的傅里叶谱：
H*(amu) (m = 1,2,…,M),(1)
假定 H 是实的,从而有
H*(amu) = H(amu),(2)
{H(amu)| m = 1,2,…,M } 构成多通道小波变换匹配滤波器,G(u) 经滤波后成为
H*(amu) G(u) (m = 1,2,…,M),(3)
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
40
xy-L2-uv 构成 x 方向的一维傅里叶变换器。在
SLM1上输入信号 g(x)，经 L2 的变换作用，在 uv
平面上形成它的傅里叶谱 G(u)．
在 uv 平面上放置第二个空间光调制器 SLM2，
它被分成 M个沿 u方向的带状区域，这些带状区域中分别显示具有不同伸缩因子 am 的基元函数 h
的傅里叶谱：
H*(amu) (m = 1,2,…,M),(1)
假定 H 是实的,从而有
H*(amu) = H(amu),(2)
{ H(amu)| m = 1,2,…,M } 构成多通道小波变换匹配滤波器,G(u) 经滤波后成为
H*(amu) G(u) (m = 1,2,…,M),(3)
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
41
uv - L3 -构成像散系统在 子午面 (vz平面 )内,L3使 uv 平面成像在上 ;
在 弧矢面 (uz平面 )内,柱面镜没有作用,uz 位于球面镜的前焦面, 位于其后焦面,构成沿 u 方向的一维傅里叶逆变换,参见图 6.8．
图 6.8(a)在子午面内构成成像系统
(b)在弧矢面内构成一维傅里叶逆变换系统第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
42
由于成像作用,在平面上沿?方向相应形成 uv 平面上各带状通道的像,对于第 m 个通道,
由于沿 u 方向的傅里叶逆变换作用,得到
(a,am) =? H*(amu) G(u) exp(i2?a u) (4)
由节 6.3(4)式,有
(5)
在图像处理系统中,将 CCD输出的信号除以
(am )1/2即得到小波变换,
在 中,伸缩因子 a 是分立的,由一组 M个滤波器 引入处理器；而位移因子
b 则是连续的,与输入平面的坐标 a成正比,
) M,2,1,m ( )a,a(a}g{W mma,a m
}g{W a,a m
),a(H m u
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
43
光学 Morlet 小波变换实例实验采用下图所示的光学多通道小波变换相关器．两个柱面透镜的焦距为 240 mm，球面透镜焦距为 210 mm．相关器由 He-Ne激光作光源．
实现一维小波变换的二维光学相关器小波函数是余弦 ——高斯 Morlet 函数它的傅里叶变换是
H(v)是一个实函数，表达式
H(0) = exp(-2?2?2v02)≠ 0，
但是随余弦函数中的频率 v0 的增加而迅速下降．
例如,当 2?v0 = 5,σ= 1时,H(0) = 3.7× 10-6，
因此小波变换的存在条件在一个非常好的近似范围内被满足．
2 te x p)t2c o s (
2
1)t(h
2
2
0
]})(2e x p [])(2{ e x p [21)(H 20222022
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
45
光学多通道小波变换相关器中光学小波变换的实验结果第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
46
6.5 用多通道匹配滤波实现二维小波变换
6.5.1 单通道小波变换系统一维小波变换定义,
二维小波变换也可类似定义：
(1)
为了简单起见，令 ax = ay = a，即 x,y方向按相同的尺度加以缩放，得到
(2)
)b(g)ab(ha1dx)x(g)a bx(ha1)}x(g{W *b,a
dd),(
a
b
,
a
b
h
aa
1
)b,b(
a
b
,
a
b
h
aa
1
)}y,x({W
y
y
x
x*
yx
yx
y
y
x
x
yx
b,b,a,a yxyx
dd),(a b,a bha1 )}y,x({W yx*b,b,a yx
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
47
在频域中,(1)式变成 匹配滤波的频域表达式 (3)
匹配滤波可以用标准的 4f系统实现,如图 6.9所示,
图 6.9 用 4f 系统实现小波坐换将二维信号函数?(x,y)经过 SLMl 输入系统,
则在 Ll的频谱面上将出现它的谱?(u,v)．
d)dbb(2ie x p [),()a,a(Haa}{W yxyx*yxb,b,a,a yxyx vuvuvuvu
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
48
在谱面上放置第二个 SLM2,将匹配滤波函数 H*(axu,ayv)通过 SLM2 对?(u,v) 进行滤波,
则形成 H*(axu,ayv)?(u,v),再经过第二个透镜
L2,在输出平面上得到它的傅里叶逆变换,
由上式可知,H*?的傅里叶逆变换,即信号?(x,
y) 的小波变换
d)dbb(2ie x p [),()a,a(Haa}{W yxyx*yxb,b,a,a yxyx vuvuvuvu
)}y,x({W yxyx b,b,a,a?
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
49
我们注意到,位移因子 ( bx,by )是与输出平面的坐标对应的变量,但伸缩因子
( ax,ay )却是给定的,亦即我们只能对给定的伸缩因子 ( ax,ay ) 实现小波变换,不同的 ( ax,ay )的变换只能通过依次输入不同的匹配滤波函数 H*( axu,ayv ) 来实现,速度很慢,发挥不了光学系统并行处理的优越性,
6.5.2 用 Dammann 光栅进行多通道相关处理由于 comb(x)与任意函数?(x)相乘,等于对?(x) 的抽样； comb(x)与?(x)的乘积再进行傅里叶变换,在谱面上得到 comb(u)与?(u)的卷积,相当于谱项?(u)在 comb(u)
中各个? 函数的位置上重复出现,只要?(x) 的带宽是有限的,则我们总可以通过足够密集的抽样手续,使各谱项在频域内互相分离,即 对于给定的输入信号?(x,y),可利用梳状函数在谱面上复制出一系列谱函数?mn(u,v)，
m,n = 0,± 1,± 2,…,如果每个?mn 代表一个处理通道,就可以实现多通道的并行处理,在这里,关键的器件是梳状函数器件,又称 Dammann光栅,
一维 Dammann光栅的透过率函数可表为
(4)
式中?是栅线宽度,d 是空间周期,见图 6.10(a)．
mdxr e c tdxc o m b*xr e c td1)x(t
m
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
51
若将 Dammann光栅放在 4f系统的输入面上,则在谱面上得到
T(u) =? sinc(?u)comb(du),(5)
它是一个被 sinc函数调制的梳状函数，见图
6.0(b),若?足够小，则 sinc(?u) 变化缓慢，
在中心附近的一些?函数幅度接近于 1．当然,? 越小，透过的光能量就越少．
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
52
把两个同样规格的一维 Dammann光栅相对旋转 90o并叠在一起,就得到二维 Dammann光栅,其透过率函数可表为
(6)
在?很小时,上式可近似表为
(7)
显见 二维 Dammann光栅是一个? 函数点阵,它们分布在间隔为 d 的栅线的交点上,这些地方的透过率为 1,其余地方的透过率为 0．
dyc o m b*yr e c tdxc o m b*xr e c td 1)y,x(t 2?
dy,dxc o m bd1)y,x(t 2?
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
53
在上图所示的 4f系统的输入平面上放置信号
(x,y),把二维 Dammann光栅紧贴在?(x,y) 上，
在光栅后的光场复振幅分布为
(x,y) comb(x/d,y/d) /d2 (8)
在 uv 平面上得到它的傅里叶谱
(9)
式中 fo = 1/d,(10)
d
1)nf,mf(
d
1
)dd(c o m b*)()(
m n
n,m2
m n
oo2
vu
vu,vu,vu,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
54
m.n =?(u - mfo,v - nfo) (11)
是一系列中心位于 Amn=(mfo,nfo) 处复现的谱项,如 设计一个匹配滤波器列阵,位于 Amn
处的复数透过率为 H*m.nexp[i2?(upm+vqn)]，
亦即：
(12)
滤波器中各谱项是按角度编码的，相当于用不同方向传播的平面波 exp [ i2?(upm+vqn)]
作为参考光制成的全息匹配滤波器，通常是计算机生成的全息匹配滤波器 (CGH)．
e)]nf(a),mf(a[H)(F
m n
)qp(2
omom
* nm vuvuvu,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
55
经滤波后,再通过 4f 系统中第二个傅里叶透镜 L2的作用,在输出平面上得到傅里叶逆变换：
式中 Cmn 是常数系数,(13)
m n
qb,pb,a,amn2
m n
nymx
n
ny
m
mx
mn2
m n
)]qb()pb([2i
nm
*
mn2
)]qb()pb([2i
oo
m n
onom
*
2yx
)}y,x({WC
d
1
)qb,pb(
a
qb
,
a
pb
hC
d
1
dde),()a,a(HC
d
1
dde)nf,mf(
)]nf(a),mf(a[H
d
1
)b,b(
nymxnm
nymx
nymx
vuvuvu
vuvu
vu
vu
vu
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
56
亦即频域中的相移形成输出平面上的位移，以 ( pm,qn )为中心形成一系列在空间相互分离的项，每一项都代表一个不同的伸缩因子 ( am,an )的小波变换，位移因子则由以 ( pm,qn )为中心的坐标来表示,
这样一来，我们就用 多通道 4f 系统实现了二维光学小波变换，它是缩放因子的分立函数，是位移因子的连续函数．
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
57
6.5.3用体全息存储器进行多通道相关处理重铬酸明胶 (DCG) 或 光 折 变 晶 体
(BaTiO3,LiNbO3,BSO等 )都可以做成体全息存储器,并可以用角度对全息图进行编码,这一性质也可以用于光学小波变换,
参见图 6,11．
将缩放后的母函数 hmn = h( x/am,y/an ) (14)
经 SLM1输入系统,则在频谱面上产生了它的傅里叶变换
Hmn = am an H ( amu,anv ),在谱平面上放置体全息记录器件,例如光折变晶体,并用倾斜的平行参考光
Rmn = exp[ -i2?(upm+vqn)] (15)
照射,与物光 am an H ( amu,anv ) 相干叠加,形成全息图,一系列 { hmn } 依次输入系统,并用不同角度的参考光 { Rmn } 编码,最后得到,(16)
式中 第三项正是我们所需的小波变换匹配滤波函数,
eHaaeHaaHaa1{
RHaa)(F
m n
)qp(2
mnnm
)qp(2*
mnnm
2*
mn
2
n
2
m
m n
2
mn
*
mnnm
nmnm
vuvu
vu,
当我们将所分析的信号?(x,y)通过 SLM1输入系统,
它的谱?(u,v) 经过 F(u,v)的滤波后得到
(17)
我们略去了与匹配滤波无关的项,再经 L2 的傅里叶逆变换,在输出平面上得到同样得到以 (pm,qn)为中心的一系列空间互相分离的小波变换．
)e()aa(Haa
m n
)qp(2
nm
*
nm
nm vuvu,v u,
m n
qb,pb,a,a
m n
nymx
n
ny
m
mx
)]qb()pb([2i
m n
nm
*
nmyx
)}y,x({W
)qb,pb(
a
qb
,
a
pb
h
dde),()a,a(Haa)b,b(
nymxnm
nymx
vuvuvu
vu
(18)
(16)式前二项经逆变换重合于原点，但第四项则是以 (-pm,-qn) 为中心的卷积项，在输出平面上与相关项关于原点对称，可称为,鬼像,
(ghost image)，因此在角度编码时要避免这些鬼像对其他小波变换项的干扰．
目前全息存储密度已做得相当高，因此在晶体中可以存储相当多个按方向编码的小波变换匹配滤波器．
eHaaeHaaHaa1{
RHaa)(F
m n
)qp(2
mnnm
)qp(2*
mnnm
2*
mn
2
n
2
m
m n
2
mn
*
mnnm
nmnm
vuvu
vu,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
61
6.6 光学小波变换匹配滤波器在图像识别中的应用图像识别或特征识别,
从大量信息或背景中检测某一特定的图形或指定的特征信息,并排斥其他图形信息 。
在一般情况下,我们只知道需要检测的图形的特征,对于其他图形的特征,我们事先可能并不知道,或知之甚少,但检测系统必须排斥这些图形,图像识别系统的这种性能称为,排他性,。
光学小波变换识别系统在这方面的性能比常规的相关识别系统更强,因此有可能利用这一效应设计成有应用价值的小波变换图像识别系统 。
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
62
6.6.1 边缘增强效应图形的重要特征之一是它的 形状或轮廓 。 为了识别某一特定的图形,往往只需认定它的轮廓,而并不需要研究它的内部细节 。 轮廓就是图形的边缘,一旦图形的边缘被清晰地勾画出来,这一图形就容易识别了 。 相对于图形整体而言,边缘显然是局部,因此我们可以期望小波变换在边缘检测中有特殊的功效 。
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
63
选择,墨西哥帽,式母函数,定义为
(1)
引入 g(x,y) = exp[-(x2 + y2)/2] (2)
则有 h(x,y) = -(1/2?)?2g(x,y) (3)
由 F {?2g(x,y)} = - 4?2(u2+ v2)G(u,v),得到
H(u,v) = - 2? (u2+ v2)G(u,v) (4)
式中 G 是 g 傅里叶变换
G(u,v) = 2? exp[ - 2? (u2+ v2)] (5)
代入 (4)得到：
H(u,v) = 4?2(u2+ v2) exp[ - 2? (u2+ v2)] (6)
2
yxe x p)]yx(2[
2
1)y,x(h 2222
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
64
图 6.12给出 a = 0.6 及 a = 0.3 的墨西哥帽小波函数及它们的傅里叶谱,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
65
函数?(x,y) 的小波变换
7)：
亦即? 的小波变换是? 和缩放后的高斯函数 g 的相关的二次导数,
),(
a
y
,
a
x
g
a2
1
-
dd),(
a
y
,
a
x
g
a2
1
-
dd),(
a
y
,
a
x
h
a
1
)y,x(
a
y
,
a
x
h
a
1
}{W
2
2
*
y,x,a
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
66
由于
(8)
伸缩因子 a 正是高斯函数的特征尺度 。 g 和
相关的结果则是平滑效应,? 中比 a 小得多的精细结构都被平滑掉 。
此外,? 和 g 的相关运算的结果再二次求导,?(x,y) 中振幅不变的区域 (常数项 )
及线性变化的区域 (一次项 )都等于零,而在振幅变化的拐点两侧不为 0。 边界正是这样的拐点,
a2
yxe x p
a
y,
a
xg
2
22
),(ay,axga2 1-}{W 2y,x,a?
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
67
(a) 边界函数?(x)和高斯函数 g(x)
(b)Cr (a,x)= g(x/a)(x) (c)小波变换函数 D (x)
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
68
图 (a)：边界函数?(x),高斯函数 g(x/a)；
图 (b)：给出它们的相关,相关运算对?(x)
起到了平滑的作用,结果噪声都被平均掉；
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
69
可以 在边界内外侧看到小波变换的一对正、负峰，
它们明确指示了边界的位置,边界内、外区域内小波变换函数都是零．其最终效果恰恰是边界突出或轮廓突出，这正是我们所期望的．
在二维小波变换的情况下，伸缩因子 a是矢量而不是标量,a = (am,an)，一般情况下 am ≠an,
从而在两个方向上的平滑效果可分别加以控制．
图 (c),给出小波变换作为位移因子的函数 (差一个常数因子 )，即
2[g(x,y)(x)
待处理坦克模型图像和小波函数平面图像及小波变换边缘增强后实验结果图 3,待处理图像小波函数图像及小波变换边缘增强后图像
( d ) I = 2 4,a = 0,0 8
( e ) I = 2 4,a = 0,1 5
( c ) I = 2 4,a = 0,0 5
(f) I = 2 4,a = 0,2 0
( a ) 待处理图像
( b ) 小波函数图像待检测彩色图像灰度图像后边缘检测
RG
B
空间边缘检测
YUV
空间边缘检测待检测彩色图像灰度图像后边缘检测
Lαβ
空间边缘检测选择域值二值化处理第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
73
6.6.2 小波变换匹配滤波图像识别课题,一般情况下仅仅勾画出图像或区域的轮廓或边界是不够的,还要求认定输入图像中是否包含要求识别的目标,
在常规光学信息处理中,我们用匹配滤波方法达到这一目的,在小波变换信息处理中,匹配滤波方案同样适用,
设 输入 要求识别的 目标图像?(x,y) 和 另一图像? (x,y),首先对它进行二维小波变换，可采用上一节中介绍的各种方法，其结果，得到了边缘增强的小波变换谱函数 Wa,x,y{?} 和?的小波变换
Wa,x,y{?}，然后再以 W{?}和 W{?} 作为输入信号,进入第二个光学相关识别系统,例如 4f 系统，
或第五章中介绍的各种实现傅里叶变换的系统。
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
74
在谱面上,W{?} 和 W{?} 的频谱函数分别由以下二式给出：
(10)
(11)
表明：在第二个相关识别系统中,当输入函数为
W{?} 和 W{?} 时,谱面上将出现?H* 和?H*。
我们 用?(x,y)所对应的小波变换匹配滤波器
F(u,v) =?*(u,v) H(au,av) (12)
进行滤波,上式中
H(au,av) = 4?2a2(u2+v2)exp[-2?2a2 (u2+v2) (13)
vuvuvu,
vuvuvu,
vu
vu
dde)a,a(H)(a}{W
dde)a,a(H)(a}{W
)yx(2i*
y,x,a
)yx(2i*
y,x,a
经滤波后得到 |?(u,v) H*(au,av) |2 (14)
及? *(u,v)?(u,v) H(au,av) H*(au,av) (15)
再经过傅里叶逆变换,分别得到 W{?} 的自相关
a-∞∞|?*(u,v) H(au,av)|2 exp[i2?(ux+vy)]dudv
= Wa,x,y{?}? Wa,x,y{?} (16)
及 W{?}和 W{?} 的互相关
a-∞∞[?*(u,v)H(au,av)][?(u,v)H*(au,av)]
exp[i2?(ux+vy)]dudv
= Wa,x,y{?}? Wa,x,y{?} (17)
其中 W{?} 的自相关给出亮斑 (自相关峰 ),成为目标图像的特征,
由于 小波变换匹配滤波方法 给出经过 边缘增强处理的二维小波变换函数的 自相关,所以有 可能获得更好的识别效果,
由于小波变换匹配滤波器是?*H,所以 H 正是频域中的窗函数,当 a增大 时通带宽度减小,通带的中心频率向低频移动,可以 有效地抑制高频噪声,但鉴别能力下降；反之,当 a 减小 时匹配滤波器将包含更广的频段,
中心频率向高频移动,将 具有更高的鉴别能力,但容易受 噪声的干扰,因为一般的白噪声是与带宽成正比的,
a = l a = 2
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
77
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
78
6.7 光学 Haar小波变换和图形边缘探测
6.7.1 Haar变换
1910年,Haar就首次提出了 Haar变换函数：
h(x)=rect[2(x–1/4)]-rect[2(x–3/4)] (1)
它较易用光学方法实现,因此 Haar小波变换是常用的小波变换之一,h(x)的傅里叶变换为
H(v) = i2?e-i?v [1-cos(?v)]/?v (2)
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
79
6.7.2 Hanr小波变换与边缘探测
Haar小波变换 是信号函数?(x)与 Haar
函数经伸缩后的 母函数 h(x/a)相关的结果,对于一个给定的伸缩因子 a,Haar小波变换的作用如下,
(1) 在小波基元函数 ha,b(x)的正,负半周内对信号进行不加权的积分,这事实上是一个平滑或平均的过程,
(2) 将正,负半周的积分值相减,以上两个作用的综合结果是在平均的意义下求差分,
或求导数,恰恰是测出了图形的边缘,
图 6.17(a)是对一个带有低频噪声的方波 进行
Haar小波变换的结果,
在方波的两个边缘呈现一对峰,极值恰恰指示了边缘的位置,
峰峰实线 为带有噪声的 方波信号 ； 虚线 为 Haar小波变换,
图 6.17(b)是对一个高频噪声的方波 进行
Haar小波变换的结果,
在方波的两个边缘呈现一对峰,极值恰恰指示了边缘的位置,
实线 为带有噪声的 方波信号 ； 虚线 为 Haar小波变换,
峰峰图 6.17(c)方波同时具有低频和高频噪声干扰,
只要伸缩因子 a 选择得当，小波变换仍然有很高的信噪比，峰很尖锐，正确的指示了边缘的所在，充分说明 Haar小波变换的抗干扰能力．
实线 为带有噪声的 方波信号 ； 虚线 为 Haar小波变换,
峰峰第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
83
图 6.18(a)是一个带有噪声的 step函数,g(x)及其 Haar小波变换,噪声是由空间频率为
v1 = 27 及 v2 = 40 的正弦干扰信号及随机本底噪声构成．尽管信噪比很底，小波变换仍然正确地指示了 step函数的边界．
图 6.18 step函数及其 H a a r
小波变换第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
84
图 6,18(b)是信号函数 g(x)和小波函数 h(x)在频域中的行为,实线表示 G(v),虚线表示 H(v)．
G(v)的 主峰主要由 step函数的直流和低频 成分贡献,而 旁瓣 则主要是 step函数 在 边界 的跃变部分贡献的,
零频、低频主峰旁瓣
H(v)频率窗恰恰开在旁瓣干扰信号的主频 v1,v2远在频率窗外第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
85
值得一提的是,两种干扰信号的主频
v1,v2远在频率窗外,对变换没有贡献,而随机本底噪声的谱是宽带的,对结果的影响也不大,
说到底,在空域中干扰信号和随机本底都是,全局,信号,小波变换的局部处理手续恰恰抑制了它们的作用,而突出了图形的局部变化 ——边界的贡献,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
86
6.7.3 二维 Haar小波变换和图形拐角测量二维 Haar变换的母函数定义 如下：
h(x,y) =rect(x - 0.5,y- 05) + rect(x + 0.5,y + 05)
+ rect(x + 0.5,y- 05)+ rect(x- 0.5,y+ 05)
二维 Haar变换母函数如图所示．用 h(x,y)
构作的二维 Haar小波变换，
在探测与 x 轴或与 y 轴平行的边缘时都为零，但测量任意一段既不和 x 轴平行又不和 y 轴平行的边缘时不为零，因而 该变换特别适用于探型图形边缘的
“拐角”,
-1 +1
+1-1
x
y
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
87
图 6.20(a) 给出一个字母,T”作为输入图形,
图 6/20(b) 是由二维 Haar变换测量该输入图形的输出,可以看到拐角测量的效果,
专门进行 x 方向或 y 方向边缘测量的二维 Haar变换,
见图 6.21．
用于测量拐角的 Haar函数常称为,角母函数,,
用于测量边缘的 Haar函数则称为,边母函数,,
根据图形的尺寸,边缘过渡区的宽度及噪声的频谱分布,
可以采用不同的伸缩因子,
以获得最佳的信息提取效果,
图 6.21二维图形边,角测量的 Haar小波变换母函数示意图,(a)角母函数,(b) y边缘母函数,(c) x 边缘母函数,
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
89
从一维和二维的边缘探测的讨论中,
我们已看到小波变换的特点,相对于图形整体 而言,边缘 ( 边和 角 )显然 只是局部,Haar小波变换对于图形内外强度不变的或缓慢变化的部分不敏感,即使对一些急剧变化的高频噪声,只要适当选取伸缩因子,Haar小波变换仍有很强的抗干扰能力,但在图形的边缘出现强度跃处,小波变换却出现很强的峰,正确指示边界的位置,可见,在边缘探测这一类问题中,小波变换比傅里叶变换更加适用,
6.7.4 用投影 ——掩模法实现 Haar小波变换本节中介绍一种由 Yang等提出来的在空域中实现
Haar小波变换的方法,投影 ——掩模法,
图 622 用投影 ——掩模法进行 2 维 Haar小波变换
o 为点光源列阵；?l 为小波函数掩模板；
2 为输入信号 SLM；?3 为探测器从点光源 P(x,y)辐射的发散球面波照亮?1上的小波函数掩模板,在 d1的距离上经放大后通过 Ll 投射到
2上,在 Ll 和 L2 间为平行光,没有放大作用,放大率即伸缩因子,a = f /( f – d1 )
因为 P(x,y)一般不在光轴上，使母函数发生平移,
平移量的计算可参见图 6.23．
光线 PQRO2在平面 上的位移量为零，它在 st 平面即?1 平面上的位移量,
bo = O1M + MQ = xd1/f +?/a = x(d1/f + d2/fa ) (5)
式中 a 为伸缩因子,MQ= NR/a =?/a,NR = xd2/f.
bo 放大 a 倍后换成平面即?2 平面上的位移量
b = boa = x(ad1/f + d2/f) = x[d1/(f -d1)+d2/f] (6)
整个光束对于?3 上一点 P’(x’,y’)的贡献为积分
(7)
由 CCD转换成电信号送入数据处理系统,在其中补上因子 a –1/2 即成为小波变换式
(8)
倘若光源是 点光源列阵,我们得到 抽样的小波变换函数 ；若光源是 均匀的连续面光源,我们得到 连续分布的小波变换函数,将掩 模板沿轴向平移,可以改变伸缩因子 a的值,
dd),(a b,a bh)b,b(' yxyx
dd),(a b,a bh
a
1 }{W yx*
b,b,a yx
第 1节第 2节第 3节第 4节第 5节第 6节目 录第 7节第 6章
2009-7-19 光学信息处理
94
6.7.5 偏振编码掩模 ——投影方案中,非相干光 学信息处理用光强大小来表示信号,只能给出 非负的分布,
这显然 违背了小波变换的相容性条件,
若使用 相干光源,可用偏振编码来实现 Haar
变换函数的双极性,
设光源辐射的线偏振光沿 45o方向振动,可以用 Jones矩阵表为
(9)
在掩模中,+1的部分用水平方向的偏振片构成,
用 Jones矩阵表为
(10)
11J? o
00 01P
-1部分用垂直方向的偏振片构成,表为
(11)
输入的光波经过 +1和 -1部分后,分别变成
10 00P
11
11
2
1
P
45
SL MyxJ
J
J
1
0
1
1
10
00
J
P
J
0
1
1
1
00
01
J
P
/4-
o
o
o
表为方向的偏振分析器，可放置沿一前入信号方向的线偏振光。在输和”分别是沿和’
”及
’
掩模板的制作。
编码技术实现了，这样一来，我们就用的相位差为方向的偏振光，它们和分别是沿和
”
及
’
变成振分析器，部分的光波，再经过偏及经过掩模中
H aar
13545J
J
J
1
1
2
1
1
0
11
11
2
1
J
P
J
1
1
2
1
0
1
11
11
2
1
J
P
1-1
oo
-
4/
4/
Wavelet transform and its use
in edge detection
Feijun Song
1 China Daheng Corporation,P.O,
Box 9671,Beijing 100086,China
Suganda Jutamulia
2 Blue Sky Research,3047 Orchard
Parkway,San Jose,CA 95134,USA
ABSTRACT
Application of the wavelet transform to edge
detection is discussed in general,To propose two new
wavelet transforms for edge detection,Modified
Haar’s wavelet transform (MHWT) and Haar-
Gaussian wavelet transform (HGWT).To analyze the
performance of MHWT and HGWT for edge
detection in both space and frequency domains.
Application for lead inspection of surface mount
devices (SMD) in electronic industry by using MHWT
and HGWT is presented in detail.
Keywords,Wavelet transform,Machine vision &
automatic inspection,Edge detection,Lead inspection,
Surface mount devices (SMD)
1,IN TRO DUC TIO N
An i m age g en er ally cont ains s om e slow l y var yin g are as and
sh ar ply var yin g r egi on s,S har p vari eties,si n gula rities o fte n carry
inter esti ng and imp o rta nt i nfo rm atio n.
EDG E D ET ECT ION (I)
Edg es of a n i ma ge a re gene rally local ized by si ngul arities an d i rreg ula r it ies.
Edg e det ecti on is one of the basi c oper ations of al l art ifi cial,com puter,robotic,
or m ach ine visio n syst em s.
Gradi ent est im ations are cl assi cal mo dels which giv e satisf act ory result s in the
case s of sli ghtl y noi sy im ages,
For Noi sy im ages a low -pass fil te r shoul d be inc lud ed,Howe ver,its perf orm anc e
is grea tl y depende nt on the smo othi ng scale fact or of the im pul se respo nse,If the
edge to be d etec te d is shar p,one can choo se a smal l im pul se re sponse,
EDG E D ET ECT IO N (II)
Bu t if th e ed ge is sm oo th,a la rge im pu lse re spo nse will gi ve m ore
sat isfa ct ory re sul ts.
Th er ef ore,sin gl e sm oo th in g scale of im pu lse re spo nse is oft en
in suff ic ie nt in m an y ap pl ica tio ns,an d m ul ti-sc al e or va ria ble im pu lse
resp on se fo r e dg e d et ec ti on is ve ry at trac ti ve.
WA VEL ET TR A NSF O RM
M all at [ 1] su g ge ste d t hat th e w a vel et tr an sf orm is very c l os ely rel a te d
to a m ult i- s cal e i mp uls e r es po ns e for e d ge d etec ti o n,
An at tr acti ve f eatu r e of t he w a vel et tra ns fo rm is th at it ca n wel l l oc ali z e
th e re gu la rit y of sign als in b oth spa c e an d f req u en cy d om ain s,By
va ryin g it s dil ati on p ar a met e r,th e wa v el et tr an sf orm a cc o rdin gl y
ch an g es i ts s mo oth in g sc ale f act or au to m ati c all y.
2,NE W W AV E LE T T RA NS FO RM S
H aa r’ s wa ve l et is co ns t ru cte d by t wo re ct an gu l ar fu nc t i on s as ex pr es se d i n
Eq,(1) an d sh ow n in Fi g,1,
4
3
2
4
1
2 xr e c txr e c txh
,(1 )
wh er e r ect ( x ) i s defi n ed a s
2/1||0
2/1||1
x
x
x
H AA R W A VE L E T
Its Fou rier tran s for m i s
)2(,
c o s1
2
i
ieH
,
It h as a ze ro dc c om po n ent,
Fi g.1,H aa r’ s wa ve let (a ) a nd it s s pe ct ru m (b )
x
h ( x )
1
( b)
H ( )?
.
( a)
TWO NEW WAVELE TS
We pr opos ed tw o new wav ete s,Mod ified Haar ’ s wavel et and Haar-Gaussian
wav el et,[4,5],d efine d as:
)0(' qs
s
qx
rec t
s
qx
rec txh
s
,(3)
)0(e x pe x p"
22
qs
s
qx
s
qx
xh
s
(4)
q,off set param e t er s,di l at i on fact or
Modifi ed H aar’s w av elet an d Haar -Gaus si an wavel et are b ot h anti-symme tri c as
show n i n Fig,2.
Wit h an ad ditional of fset
para met er q,m odifie d Haa r’ s
wav el et and H aar-Gau ssian
wav el et ar e ver y powe rful in
noisy edg e dete ction and
featur e ext raction,
Fig,2,Modifi e d Haar ’s wave le t
( sol id line ) and
Haar -Gaus sia n wave let
( das hed line)
- q q
3,NOISY EDGE DE TECTIO N
Math em at ic al mo de l of edge functi on is expr essed i n Eq,(5)
2
1t a n
1
x
xf
,(5)
Wi dth of ed ge? is
given by
,( 6)
The ed ge functi on is
depict ed in Fi g.3.
Fig,3,Model of an edg e
with w id th?,
NOISY EDGE FUN CTION (I)
No is y edge functi on is furt her s i m ulated as f ollow s.
r n dxx
x
xf
2211
1
2s i n2s i n
2
1t a n
,(7)
Last t hree ter ms repr esent the sinu soid int erferen ce and rando m n ois e,
Si mu lat ed nois y edge functi on is depict ed i n Fi gs,4(a),(b),(c) an d ( d),t ogether with
the H aar-Gaussi an wavelets of various d il ations s.
The vertical li nes ( da shed) specif y t he sp at i al win do w s of the wavelet,
( a) ( b)
(a) s=0.63 (b) s=0.42
NOISY E DGE FUNCTION (II)
Fig,4,No is y e dge functi on a nd Haar-Gaus sian w avelets with various dil at ions,
(a) s = 0,63,2s = 9.5? ; ( b) s = 0,42,2s = 6,4? ;
(c) s = 0,10,2s = 1.6? ; (d ) s = 0.053,2s = 0,8?
2s,w idt h of th e spa ti al window, is de fi ned as follow s.
2/1
2
","
","
2
xhxh
xhxxh
s?
,(8)
where ( f,g ) repres en t s the inn er product o f f un cti on s f and g,
(c) s=0.01
(d) s=0.05 3
HAAR-GAU SSIAN WA VELE T TRA NSFORM OF TH E
NOISY EDGE FUN CTION (I)
Ha ar- G aussi an wa vel et tr ansf or m ( HGW T ) of t he noi s y edg e
functi on f ( x ) is a cor r elati on of t he sign al wit h the wa vel et:
d
s
x
hfxhxfxfW
ss?
)("
*
,
(9)
HGW T for diff er ent di lat i ons s (wh i c h is equi val ent t o th e
smo oth in g sc ale of im pul s e r es pons e dis c uss ed pr e vio usl y) ar e
depi ct e d in F ig,5.
M axi ma of th e W T in di cat e
th e locati on of the e dge,
Fi g,5,HG WT of the no i sy edge
fo r di ff erent dil a ti on s,
(a) s = 0.63 (cu rv e wit h t he
br oad est peak),
(b) s = 0.42 (cu rv e wit h t he l ess
br oad peak),
(c) s = 0.10 (cu rv e wit h nar ro w
peak and sm al l ri pp l e s),and
(d) s = 0.05 3 (cu rv e wit h sev ere
ri pp l es),
HAA R-GA USSIAN WAVELET TRANS FORM OF THE
NOI SY ED GE FUNCTI ON (II)
Charac t eris tic d i me nsi o n o f the noise is gi ven as:
2.0/1
,(10)
,th e freq ue nc y o f the sinu soi d al n o ise.
,
s,d ilat io n of WT
2s,wid th of sp atial wi n d o ws
,char acteriz ed n o is e di mens io n,
When th e wid t h o f the spa tial w i n do w is clo se t o th e ch aracteristic
d imension of t he n o ise multiple ma xim a a re found in the wavelet
transfor ms and the fal se response w o ul d oc cur.
4,FREQUEN CY DOMAIN ANALYSI S
The edge detection based on the wa velet transf orm in space domain can
be analyzed in frequency domain as we ll,The frequency spectra of both
the noisy edge function and the Ha ar-Gaussian wa velet of differ ent
dilations s are depicted in F ig.6 (a) to (d).
S 2s
0.63 3.0
0.42 2.0
0.10 0.50
0.053 0.25
0.2
Tab le 1
Fig.6,Spectra of noisy edge functi on (sol id line),H aar-G aussian
wavelet (dashed li ne) and frequency fi lt ering window (Verti cal li nes),
(c) s = 0.12
(a) s = 0.63
(d) s = 0.053
(b) s =0.42
FREQUENCY DOM A IN ANAL YSIS ( I)
Ha ar- Gaussi a n
wa vel et acts a s a
low fr e quen c y
ban d- pass fi lt er
wi th no dc
comp onent i n
fre que ncy dom ain,
The s pect r um of the r a mp edg e fu nct io n consist s of a hi gh dc
comp on ent a nd a se ri e s of low f re quen cy l obe s,The s harp est vari ati on
of t he il lu mi nant cur ve (t he e dg e) c ontri but es mo stl y t o th e l ow
fre que ncy l obes,
The ba nd-pa ss fi lt er of th e Haa r -Gau s s ia n wav e let block s th e dc
comp on ent but pas s es th e low fr equ enc y l obes belo ng t o th e e dge
sign al a nd th us e nh an ces and ext ract s th e e dge tr ansi ti on
in form ati on.
FR EQ UEN CY DO M A IN ANAL Y SIS (II )
Q q
FR EQ UEN CY DO M A IN ANAL Y SIS (II I)
In Fi gs,6( a) a nd (b),th e ban d-pas s fil t er of th e HGW T blocks th e tw o hi gh-
fre quen cy pea ks of th e sin usoi d noi se fun cti ons,These tw o pe aks a re fa r
bey ond t he cut -of f f r equ ency of t he ban d-pass fil ter,
Al so,t he spe ctru m of ra ndom noi se is broa dl y spr eadin g,havi ng le ss in flu enc e
to t he wave let tr a nsf or m.
(a) s = 0.63
(b) s=0,42
In co nt ra st t o Fi g.6 (a ) an d (b ),t he fr eq ue nc y
wi nd ow of Fi g.6 (d ) co nt ai ns t he no i se
ch ar ac t er i st i c fr eq ue nc y,whi ch pr od uc es t he
m ul ti- m ax im a i n t he sp ac e do m ai n an d m i gh t
re su l t i n som e fa ul t l oc ati ons i n ed ge
de t ec t i on,
Fi g.6 (c ) re pr es en t s a cr i tic al
si t uati on whe re t he pe ak s j us t l oc at e at
t he m ar gin of t he fr eq ue nc y wi nd ow,
Th ei r co nt ri but i on i s i nsuf fi c i en t t o
i m pa ct th e ed ge det ec tio n.
Fi g,6 (c)
Fi g 6 (c )
Fi g,6 (d)
FR EQ UN CY DO M AIN AN AL Y SIS ( IV)
CO NC L UT IO N
The fr eque n c y dom ai n a nal ysis i s i n a g o od a gre e m ent wi th t h e spa c e d om ain a n alysi s.
In t he fina l a nal ysis th e r an d o m an d sin u soid al n oise h av e in fin it e e xt en sion s in sp ace,
but th e wa v el et tr an sfo rm loc ali z es th e si gn al in b oth s pa c e an d fr equ en cy dom a i n s,It
is th ere for e i ns e nsit iv e t o the n ois e a nd abl e t o accur a tel y in di cat e th e lo c ati on of th e
edg e.
5,AP PL I CA T IO N T O LE AD INS P EC T IO N O F
SU RF AC E M O UN T D EV I CE S
BA CK G RO UND
In t he last d e ca de th e v olu me p ro d ucti o n of s urfa c e mo u nt de vic e ( S M D ) has gr o w n
trem e nd ou sly a n d t h er e is a t re n d t o r epl ac e t h e th ro u gh -hole c o m po n e nts b y S M D.
T HR OU GH- HOL E C OM P ON E N T S A N D
SUR F A C E M OU N T D E V I C E ( SM D )
S u r fac e m ou n t d e v i c e (S M D ) T h r o u gh -h o l e d e v i c e
Th i s ne w t e c hn ol og y p r od uc e s b e ne f i t s i n
M i n i at u r i z at i o n
M u c h l ar ge r i n t e gr at i o n
c os t r e d u c t i o n,a nd
d e s i g n f l e xi b i l i t y,
COP L ANA R I T Y,P I T CH AND F OO T ANGL E OF S M D
H ow e ve r,t he f ol l ow i ng f a c t o r s of S M D h a ve put m uc h s t r e ss on t he p r oduc t i on
t e c hnol ogy,v e ry f i n e r pi t c h ( <0.5 m m ) a nd m u c h h i g h e r pac ka ge l e ad c o u n t s
( m o r e t h a n s o m e hund r e ds l e a d s pe r c h i p ),
To a s s u r e r e l i a bl e s ol de r c ont a c t s,i t i s r e q u i r e d t h a t t he l e a d c opl a n a ri t y a nd
pi t c h s houl d b e 1 00 m i c r ons,The f oot a n g l e s houl d b e l e s s t h a n 10 de g.
The I C m a nuf a c t u r e r s a r e r e s pons i bl e f o r a c c u r a t e l y i n s pe c t i ng t he 3- D
pos i t i o n s of t h e l e a ds f o r a l l ou t goi ng c o m p one nt s,
C
P
C o p l an ar i ty ( C ),p i t c h ( P ) an d fo o t an gl e (? ) o f S M D
DUAL SH A DO W TE CH NIQ UE F O R LE AD INSP ECT IO N
S everal te chni qu es h ave b een pr o pos ed f or th e le ad i ns pect io n o n th e pr o du ctio n
line,T h e du al sh ad o w te chni qu e [6 ] is cos t-effec t i ve,r elia ble,an d acc ept ed b y
most m anuf act ur ers,
A s che matic di agr am of th e lead i n s pect ion syst em is sho w n in Fi g,7.
An S M D is p ut on th e ped estal,All its l ead s do n ot t ou ch th e to p pl an e of a
defu sin g scr een s o as t o m ake a n on des tru ctiv e t estin g sy ste m,
The col lim ate d beam fir st illu min ate s the l ead s o f on e sid e from l ef t,and then
the se cond col limat ed beam illumi n ates th e sa me lead s from r i gh t,
Tw o sh ad o ws ar e pr o je cte d on t he su rfa ce of a de fu sin g scr een.
The s hado ws ar e th en im ag ed by a le ns o nto a CCD an d grabbe d and
pr oces s ed b y a d ata pro ces sin g sys tem,
L E A D I N SP E C T I ON SYST E M
F i g,7,S c h e m at i c d i ag r am o f th e l e ad i n s p e c t i o n o f S M D
b as e d o n d u al s h ad o w te c h n i q u e,
c ol l i m a t e d be a m s c ol l i m a t e d be a m s
l e f t pr o j e c t e d i m a ge
d e f u s i ng
ed
s c r e e n
CCD
l e ns
p e de s t a l
s c a nn i n g
l i ne of C C D
r i g h t pr o j e c t e d i m a g e
B ASI C M E ASUR E M E NT OF L E AD I NSP E CT I ON
A b as i c m e as u re m e n t of t h e i n s pe c t i on is t o ac c u r at e l y de t e rm i n e t h e
pos i t i o n s of t h e c e n t e rs of l e ads i n a s c a nn i ng l i ne of C C D,
A n e d ge de t e c t i o n t e c hni q ue i s g e ne r a l l y a ppl i e d t o de t e r m i ne t he f ro n t a l
an d t ra i l i n g e dge s of t he s h a dow of a l e a d,The c e n t e r po si t i o n c a n t he n b e
de t e r m i ne d b y t he pos i t i on s of e dge s.
NUMERI CAL DIF ERENT IAT ION A ND W AVELE T TRA NSFORM
Numerical dif ferentiation is often unab le to use for the following reasons:
1,The edge is a ramp that canno t be d efi ned definitel y,
2,T he detect ed gray scale c urve o f leads is corrupted by nois e,
3,T he cu r ve flu ctuates du e t o t he aging of the li gh t source an d especiall y du e t o t he
va ri ation of the backgro und brig htness.
To overcom e t he pro blem,w e ap ply the newly prop osed MHWT and H GWT,T he
image of a lead,wh i ch likes a va l l ey in the gray s cale cur ve (ima ge s of le a ds )
shown i n Fig,8,is al m ost s ymmetric,
Therefore w e do not de t ec t the edges of t he valleys.
Ins t ead,w e de t ermine t he symm etri ca l center s of the valleys dire ct l y using the
MHWT and HG WT,The algorithm is called sym me t ri c evaluati on,
s
Fig,(8)
O FF SE T A ND DIL AT IO N P AR AM ET ERS
In t he WT,th e o ff set p ar ame te r q re p r es e nts t he offs e t of th e e d ge t o t h e
sym m etri c al cent e r,a nd w e s ele ct t he dil a ti on pa ra m eter s to ma t c h th e wi dt h of
the ra m p o f edg e ap pr oxi mat ely,
q q 2s
AD AP TI VE L YL Y S M OO TH I NG S CA LI NG
Th e wa ve let tra ns fo rm pr od uc es m ul ti-s cale im pu lse res po ns e.
W e ha ve fu rth erm or e de sig ne d a sp ecia l m eri t fu nc tio n to
ev alu ate th e in sp ecti on pr oc ess in g,W he n th e wi dt h of th e
ed ge is va rie d du e to so m e re aso ns,th e di la tio n pa ra m eter s is
ch an ge d to m at ch th e wid th of th e ed ge ad ap tive ly,ca lled
ad ap tive ly sm oo th in g s ca lin g.
EXP ERIM ENT RES ULT
Fi g,8,The de te ct ed gr ay sc al e cu rv e of th e l ea d sha do w (l i ne
wi th + si gn s) an d i ts m od i fie d H aa r wav el et tr an sf or m
(dashed line) and H aa r- G au ssi an wav el et tr an sf or m (s oli d
l i ne ).
Bec a us e b oth mo difi ed Ha ar ’s an d Haar - Gau ss ian
wav el ets a r e an ti -sy m m etri c,a n d t h e vall e y of th e si gn al is
alm ost sym m etri c t h e wa v elet tr a nsf or ms of t h e gra y sc ale
curv e ( c orr elati on o f si gn al a n d th e w a vel et ) pro d u ce
dire ctl y th e z er o-c ro ssi n gs at t he c ente rs o f th e v a ll e ys
(i ma g e of l eads ),
T he m a xim a an d mi n im a o f th e w avel et t ran sfor ms
re pre sen t th e e d ges of th e lea d sh ad ows,an d t he zero-
cr ossi n g s accu r atel y in di cate th e c en ter s of l eads,
A P P L I C A T I O N O F W A V E L E T T R A N S F O R M T O
L E A D I N S P E C T O N O F S M D
T he w a ve l e t t r a ns f or m s a r e r e p r e s e nt e d b y
T he w a ve l e t t r a ns f or m i s a c o n t i n u o u s f u n c t i o n o f t h e s p a c e v a r i a b l e
x,a nd i s a d i s c r e t e f u n c t i o n o f t h e f r e q u e n c y v a r i a b l e ( d i l a t i o n p a r a m e t e r ) s,
A s p e c i a l m e r i t f u n c t i o n i s d e s i n e d t o o p t i m i z e t h e d i l a t i o n p a r a m e t e r s
a s di s c u s s e d b e f or e,
d
s
x
hfxhxfxfW
ss?
)("
*
n
ssss,...,,
21
P H O T O G R A P H Y O F C O N S T R U C T E D S Y S T E M F O R L E A D
IN S P E C T IO N
T h e p h o to g ra p h o f t h e r e al se tu p i s s h o w n i n F ig,9,
T h e i n st ru m e nts h a v e b e e n p ut in to m a n u fa ctu ri n g,
CO NCL U DING RE M A RK S
W e ha ve pro po sed ne w m od ifie d Ha ar ’s wa ve le t tra nsfo rm an d
Ha ar -Ga ussi an wa ve le t t ra nsfo rm,
Th e pro po sed wa ve le t tra nsfo rm s ca n be ap pli ed to th e ed ge
de te cti on,
Bo th wa ve le t tra nsfo rm s ar e al so pro ve d ve ry po we rfu l in find in g
th e sym m etri c ce ntra l po siti on of no isy ba rs,whi ch can be ap pli ed
to th e r eal le ad insp ec tio n o f su rfa ce m ou nt devi ce s (SM D).
Th e in stru m ents ha ve be en put in to m an ufa ctur in g an d ap plie d in
th e p rod uc ti on lin es.
DIF FER ENCES BE TWEEN FOUR IER TRAN SF ORM AN D W A VEL ET TRAN SFOR M
T he Fouri e r transform (FT) e xtrac ts the p e riod ic signals on infinite duration,while t he
wa velet transform (WT ) look s a t the window specified by its dilation para m e t e r s a nd
c a n be shift e d by its transla tion para m e ter x,T h e FT i s suitab le for proce ssing
signals on infinite or long r ange duration,but the WT is c a p a ble of proce ssing
transition o r short -term signals.
For a c tual a ppli cation the a ppare nt non- z e ro c o m ponents of a wa v e let only e xist in a
very sm a ll region in both spac e a nd freque ncy dom a ins,An a ttrac tive feature of the
WT is that it lo c a liz e s the proce ssing of signal in both spac e a nd freque ncy dom a ins.
T here fore WT is a space - and freque ncy localizati on operat or,
DIFFE REN CES BET WEE N FO UR IER TRA NS FO RM AN D W A VELET
TRA NS FOR M
The fr eq u e n c y v aria bl e? a n d s p ace vari a ble x ap pe a r i n th e F T a nd
inv erse F T,re sp ec ti ve l y
.)2e x p ()()(
,)2e x p ()()(
dxjFxf
dxxjxfF
The WT p r o vid es in si gh t th at simu lt ane ousl y com bin e s th e
vari a ble s in b oth sp a c e and fr eq u ency dom ain s,i.e,th e
tran slati on p aram eter x a n d d il atio n p aram eter s ap pe a r i n th e WT
simult an e o u sly,
.)()}({
*
d
s
x
hfxfW
s
