Skip to content

Latest commit

 

History

History
1707 lines (1268 loc) · 79.9 KB

4 图像增强.md

File metadata and controls

1707 lines (1268 loc) · 79.9 KB

[TOC]

第四章 图像增强

4.1 引言

通过一些运算来突出图像的某些信息,如轮廓等.但通过运算的图像不一定逼近原图像. 图像增强:

  • 空间域法

$$ g(x,y)=h(x,y)*f(x,y) $$

  • 频域法

$$ \stackrel{F(x, y)}{\longrightarrow} \boxed{正交变换} \stackrel{F(u, v)}{\longrightarrow} \boxed{修正H(u,v)} \stackrel{G(u, v)}{\longrightarrow} \boxed{逆变换} \stackrel{G(x, y)}{\longrightarrow} $$

4.2 直方图

4.2.1 直方图的基本概念

如果将图像中像素亮度(灰度级别)看成是一个随机变量,则其分布情况就反映了图像的统计特性,这可用Probability Density Function (PDF) 来刻画和描述,表现为灰度直方图 (Histogram)。灰度直方图是灰度级的函数,它表示图像中具有某种灰度级的像素的个数,反映了图像中每种灰度出现的频率,如图4-1所示。灰度直方图的横坐标是灰度级,纵坐标是该灰度级出现的频度,它是图像最基本的统计特征。 $$ \begin{aligned} \begin{array}{|l|l|l|l|l|l|} \hline 1 & 2 & 3 & 4 & 5 & 6 \ \hline 6 & 4 & 3 & 2 & 2 & 1 \ \hline 1 & 6 & 6 & 4 & 6 & 6 \ \hline 3 & 4 & 5 & 6 & 6 & 6 \ \hline 1 & 4 & 6 & 6 & 2 & 3 \ \hline 1 & 3 & 6 & 4 & 6 & 6 \ \hline \end{array} &\Rightarrow& \begin{array}{|c|c|c|c|c|c|} \hline 1 & 2 & 3 & 4 & 5 & 6 \ \hline 5 & 4 & 5 & 6 & 2 & 14 \ \hline \end{array}\ \ 原图&&灰度直方图 \end{aligned} $$ image-20230423190956673

图4-1 灰度直方图


image-20230423191147606image-20230423191158008

灰度直方图有以下几个重要概念:

  1. Pixels:像素点数。表示图像中的总像素数,是直方图分析的基础。
  2. Level:灰度级数。表示直方图划分的灰度范围被等分为几级。Level越高,直方图越详细。
  3. Mean:平均灰度。表示图像中所有像素的平均灰度值。
  4. Count:计数。表示每个灰度级中的像素个数。Count越大,该灰度级所占比重越大。
  5. MSE(Mean Square Error):均方误差。表示直方图中各灰度级与图像平均灰度的误差的平方和,用于评估直方图的均匀性。MSE越小,直方图越均匀。
  6. Percentile:百分位数。表示小于等于该灰度级的像素占总像素的百分比。如25%的像素集中在20灰度或以下,50%的像素集中在40灰度或以下等。
  • Pixels表示图像的像素总数,是直方图产生的基础。
  • Level决定了直方图的细节程度。
  • Mean表示图像的平均亮度。
  • Count反映每个灰度级的比重。
  • MSE反映直方图的均匀性,用于评价直方图的质量。
  • Percentile显示灰度级以下的像素占比,用于判断图像的亮度分布。

设$r$代表图像中像素灰度级。 作归一化处理后,$r$将被限定在$[0, 1]$之内。在灰度级中,$r=0$代表黑,$r=1$代表白。 对于一幅给定的图像来说,每一个像素取得$[0, 1]$区间内的灰度级是随机的,也就是说$r$是一个随机变量。假定对每一瞬间,它们是连续的随机变量,那么就可以用概率密度函数$p_r (r)$来表示原始图像的灰度分布。如果用直角坐标系的横轴代表灰度级$r$,用纵轴代表灰度级的概率密度函数$p_r (r)$,这样就可以针对一幅图像在这个坐标系中作出一条曲线来。这条曲线在概率论中就是概率密度曲线,如图4-2所示。

将归一化后的图像灰度值$r$看作随机变量的主要理由有:1. $r$可以取$[0, 1]$区间内的任意值。虽然图像是确定的,但每个像素的灰度值在某一瞬间可以视为随机选取$[0, 1]$区间内的某个值。2. $r$的取值受很多随机因素的影响。如图像光照条件的变化、传感器噪声的影响、失真等会引起$r$的值产生随机变化。3. 在概率论意义上,随机变量是指数值变化不完全遵循确定规律,结果在一定程度上具有不确定性的变量。$r$满足这一定义,可以看作是概率空间内的一维随机变量。4. 假定$r$在某一时刻是连续的,这意味着$r$可以取$[0, 1]$内的任意实数值,增加了$r$的随机性。如果$r$是离散的,只能取有限个离散值,其随机性会相对较小。

image-20230423191238368

图4-2 图像灰度分布概率密度函数

下图为灰度值的曲线(仅供了解)

image-20230423191307556

$$ P(r)=\int_0^r p(r) d r \ p(r)=\frac{d P(r)}{d r} $$ 从图像灰度级的分布可以看出一幅图像的灰度分布特性。例如,从图4-2中的(a)和(b)两个灰度分布概率密度函数中可以看出: (a)的大多数像素灰度值取在较暗的区域,所以这幅图像肯定较暗,一般在摄影过程中曝光过强就会造成这种结果;(b)图像的像素灰度值集中在亮区,因此,图像(b)将偏亮,一般在摄影中曝光太弱将导致这种结果。当然,从两幅图像的灰度分布来看图像的质量均不理想。

4.2.2 直方图的性质

(1) 直方图是一幅图像中各像素灰度值出现次数(或频数) 的统计结果,它只反映该图像中不同灰度值出现的次数(或频数),而未反映某一灰度值像素所在位置。也就是说,它只包含了该图像中某一灰度值的像素出现的概率,而丢失了其所在位置的信息。

(2)任一幅图像,都能惟一地确定出一幅与它对应的直方图,但不同的图像,可能有相同的直方图。也就是说,图像与直方图之间是多对一的映射关系。如图4-3就是一个不同图像具有相同直方图的例子。

image-20230423191438660

图4-3 图像与直方图间的多对一关系

(3) 由于直方图是对具有相同灰度值的像素统计得到的,因此,一幅图像各子区的直方图之和就等于该图像全图的直方图,如图4-4所示。

image-20230423191508569

图4-4 直方图的分解

可以产生边界线

4.2.3 直方图的计算

考题:灰度直方图的计算公式

灰度直方图的计算非常简单,依据定义,在离散形式下,用 $r_k$代表离散灰度级,用$p_r (r_k )$代表pr (r),并且有下式成立: $$ p_r\left(r_k\right)=\frac{n_k}{n} \quad 0 \leq r_k \leq 1 \quad k=0,1,2, \cdots, l-1 \tag{4-1} $$ 式中:$n_k$为图像中出现$r_k$级灰度的像素数,$n$是图像像素总数,而$n_k /n$即为频数。在直角坐标系中做出$r_k$与$p_r(r_k)$的关系图形,即称为该图像的直方图。设若图像具有$L$级灰度(通常$L=256$,即8位灰度级),则大小为$M×N$的灰度图像$f(x, y)$的灰度直方图pBuffer$[0…L-1]$可用如下算法得到:

(1) 初始化

pBuffer[k]=0 ; k=0, …, L-1

(2) 统计

pBuffer[f(x, y)]++; x, y =0, …, M-1, 0, …, N-1

(3) 归一化

pBuffer[f(x, y)]/=M*N

其中,直方图的归一化是一个可选项,若不需要特殊处理可以不进行此项操作。 图4-5和图4-6分别是Lena图像、 钟楼图像与其对应的直方图。

image-20230423191605444

图4-5 Lena图像及直方图 (a) Lena图像;(b) Lena图像的直方图

image-20230423191626174

图4-6

4.2.4 直方图的拉伸

将感兴趣的部分增强

如上所述,一幅给定图像的灰度级分布在$0\leq r\leq 1$范围内。可以对$[0, 1]$区间内的任一个$r$值进行如下变换: $$ s=T(r)\tag{4-2} $$ 也就是说,通过上述变换,每个原始图像的像素灰度值$r$都对应产生一个$s$值。变换函数$T(r)$应满足下列条件:

  • $0 \leq r \leq 1$ 区间内, $T(r)$ 值单调增加;
  • 对于 $0 \leq r \leq 1$, 有 $0 \leq T(r) \leq 1$

这里的第一个条件保证了图像的灰度级从白到黑的次序不变,第二个条件则保证了映射变换后的像素灰度值在容许的范围内。满足这两个条件的变换函数的例子如图4-7所示。

从$s$到$r$的反变换可用下式表示: $$ r=T^{-1}(s)\tag{4-3} $$ 由概率论理论可知, 如果已知随机变量 $\xi$ 的概率密度函数为 $p_r(r)$, 而随机变量 $\eta$$\xi$ 的函数, 即 $\eta=T(\xi), \eta$ 的概率密度为 $p_s(s)$, 所以可由 $p_r(r)$ 求出 $p_s(s)$

因为 $s=T(r)$ 是单调增加的, 由数学分析可知, 它的反函数 $r=T^{-1}(s)$ 也是单调函数。在这种情况下, 如图4-8所示, $\eta&lt;s$ 且仅当 $\xi&lt;r$ 时发生, 所以可以求得随机变量 $\eta$ 的分布函数为 $$ F_\eta(s)=p(\eta<s)=p[\xi<r]=\int_{-\infty}^r p_r(x) d x\tag{4-4} $$ 对式(4-4)两边求导, 即可得到随机变量 $\eta$ 的分布密度函数 $p_{s}(s)$ 为 $$ p_s(s)=p_r(r) \cdot \frac{d}{d s}\left[T^{-1}(s)\right]=\left[p_r(r) \cdot \frac{d r}{d s}\right]_{r=T^{-1}(s)}\tag{4-5} $$ 通过变换函数$T(r)$可以控制图像灰度级的概率密度函数,从而改变图像的灰度层次。这就是直方图修改技术的理论基础。

image-20230423192100258

图4-7 灰度拉伸变换函数

image-20230423192108271

图4-8 $ r$和$s $的变换函数关系

在选择了变换函数$T(r)$之后,我们需要根据公式(4-5)计算新图像的灰度分布$p_s(s)$。这个公式的意思是,新图像的灰度分布$p_s(s)$可以由原图像的灰度分布$p_r(r)$和$T(r)$的求导$d r/ d s$计算得来。 也就是说,对任意一个新灰度$s$,我们首先利用$s=T(r)$求得对应的原灰度$r$。然后利用这个$r$在原图像的灰度分布$p_r(r)$中查找到其概率密度。最后,乘以$T(r)$在$r$点的斜率$d r/d s$,就得到新的灰度$s$的概率密度$p_s(s)$。 重复这个过程,就可以计算出新图像在整个[0, 1]区间内的灰度分布$p_s(s)$。 有了$p_s(s)$之后,我们就可以生成新图像了。对新图像的每一个像素,我们可以从$p_s(s)$中随机采样得到其灰度值$s$,利用$s=T(r)$再求得其对应的原灰度$r$。这样每个像素的新灰度$s$就由原来的$r$通过变换函数$T(r)$映射得到了。 最后一步是重新量化。因为数字图像的灰度是离散的,通常为0-255或0- 65535,所以我们需要将[0, 1]区间内的连续灰度$s$映射到这些离散灰度上。这个过程就是量化,通常采用线性映射,如: $$ s\longmapsto \lfloor 255\cdot s\rfloor $$ 经过上述过程,我们就得到了直方图被拉伸后的新图像。这个新图像具有扩展的动态范围和增强的对比度,能够更好地显示图像中的细节与特征。

4.2.5 直方图均衡

(重点,又称均衡拉伸算法)

直方图均衡化处理是以累积分布函数变换法为基础的直方图修正法。假定变换函数为 $$ s=T(r)=\int_0^r p_r(\omega) d \omega\tag{4-6} $$ 式中: $\omega$ 是积分变量, 而 $\int_0^r p_r(\omega) d \omega$ 就是 $r$ 的累积分布函数。 这里, 累积分布函数是 $r$ 的函数, 并且单调地从 0 增加到 1 , 所以这个变换函数满足关于 $T(r)$$0 \leq r \leq 1$ 内单值单调增加。在 $0 \leq r \leq 1$ 内有 $0 \leq T(r) \leq 1$ 的两个条件。

对式中的 $r$ 求导, 则

$$ \frac{d s}{d r}=p_r(r)\tag{4-7} $$

再把结果代入(4-5)式有

$$ p_s(s) =\left[p_r(r) \cdot \frac{d r}{d s}\right]{r=T^{-1}(s)}=\left[p_r(r) \cdot \frac{1}{d s / d r}\right]{r=T^{-1}(s)} =\left[p_r(r) \cdot \frac{1}{p_r(r)}\right]=1 \tag{4-8} $$ 由上面的推导可见,在变换后的变量s的定义域内的概率密度是均匀分布的。因此,用r的累积分布函数作为变换函数,可产生一幅灰度级分布具有均匀概率密度的图像。其结果扩展了像素取值的动态范围。

例如,在图4-9中,(a)是原始图像的概率密度函数。从图中可知,该图像的灰度集中在较暗的区域,是一幅曝光过强的照片。由图(a)可知,原始图像的概率密度函数为 $$ p_r(r)= \begin{cases}-2 r+2 & 0 \leq r \leq 1 \ 0 & \text { 其他 }\end{cases} $$ 用累积分布函数原理求变换函数 $$ s=T(r)=\int_0^r p_r(\omega) d \omega=\int_0^r(-2 \omega+2) d \omega=-r^2+2 r $$ 变换后的s值与r值的关系为 $$ s=-r^2+2 r=T(r) $$ 按照这样的关系变换,就可以得到一幅改善质量的新图像。这幅图像的灰度层次将不再是呈现较暗色调的图像,而是一幅灰度层次较为适中,比原始图像清晰,明快得多的图像。可以证明,变换后的灰度及概率密度是均匀分布的。图4-9(b)和(c)分别为变换函数和变换后的均匀的概率密度函数。

image-20230423192245548

图4-9 直方图变换法

上述方法是以连续随机变量为基础进行讨论的。当灰度级是离散值时,可用频数近似代替概率值,即 $$ p_r\left(r_k\right)=\frac{n_k}{n} \quad 0 \leq r_k \leq 1 \quad k=0,1, \Lambda, l-1\tag{4-9} $$ 式中: $l$ 是灰度级的总数目, $p_r\left(r_k\right)$ 是取第 $k$ 级灰度值的概率, $n_k$ 是图像中出现第 $k$ 级灰度的次数, $n$ 是图像中像素总数。

式(4-6)的离散形式可由式(4-10)表示: $$ s_k=T\left(r_k\right)=\sum_{j=0}^k \frac{n_j}{n}=\sum_{j=0}^k p_r\left(r_j\right) \quad 0 \leq r_j \leq 1 \quad k=0,1, \cdots, l-1\tag{4-10} $$ 其反变换式为 $$ r_k=T^{-1}\left(s_k\right)\tag{4-11} $$ $$ \begin{array}{|c|c|c|} \hline r_k & n_k & p_r\left(r_k\right)=\frac{n_k}{n} \ \hline r_0=0 & 790 & 0.19 \ \hline r_1=1 / 7 & 1023 & 0.25 \ \hline r_2=2 / 7 & 850 & 0.21 \ \hline r_3=3 / 7 & 656 & 0.16 \ \hline r_4=4 / 7 & 329 & 0.08 \ \hline r_5=5 / 7 & 245 & 0.06 \ \hline r_8=6 / 7 & 122 & 0.03 \ \hline r_7=1 & 81 & 0.02 \ \hline \end{array} $$

表4-1 64×64大小的图像灰度分布表

image-20230525141018795

图4-10 直方图均衡化处理

处理过程如下: 由式(4-10)可得到变换函数 $$ \begin{aligned} & s_0=T\left(r_0\right)=\sum_{j=0}^0 P_r\left(r_j\right)=P_r\left(r_0\right)=0.19 \ & s_1=T\left(r_1\right)=\sum_{j=0}^1 P_r\left(r_j\right)=P_r\left(r_0\right)+P_r\left(r_1\right)=0.44 \ & S_2=T\left(r_2\right)=\sum_{j=0}^2 P_r\left(r_j\right)=P_r\left(r_0\right)+P_r\left(r_1\right)+P_r\left(r_2\right)=0.65 \ & S_3=T\left(r_3\right)=\sum^3 P_r\left(r_j\right)=P_r\left(r_0\right)+P_r\left(r_1\right)+P_r\left(r_2\right)+P_r\left(r_3\right)=0.81 \end{aligned} $$ 依此类推:$S_4=0.89, \quad s_5=0.95, \quad ~s_6=0.98, \quad ~s_7=1.0$ 。

变换函数如图4-10(b)所示。

这里只对图像取8个等间隔的灰度级,变换后的值也只能选择最靠近的一个灰度级的值。因此,对上述计算值加以修正: $$ s_0 \approx \frac{1}{7}& s_1 \approx \frac{3}{7}& s_2 \approx \frac{5}{7}& s_3 \approx \frac{6}{7} \ s_4 \approx \frac{6}{7}& s_5 \approx 1& s_6 \approx 1& s_7 \approx 1 $$ 由上述数值可见,新图像将只有5个不同的灰度级别,可以重新定义如下一组符号 $$ s_0^{\prime}=\frac{1}{7}, s_1^{\prime}=\frac{3}{7}, s_2^{\prime}=\frac{5}{7}, s_3^{\prime}=\frac{6}{7}, s_4^{\prime}=\frac{6}{7} $$ 因为 $r_0=0$, 经变换得 $s_0=1 / 7$, 所以有 790 个像素取 $s_0$ 这个灰度值。 $r_1$ 映射到 $s_1=3 / 7$, 所以有 1023 个像素取 $s_1=3 / 7$ 这一灰度值。依次类推, 有 850 个像素取 $s_2=5 / 7$ 这个灰度值。但是, 因为 $r_3$$r_4$ 均映射到 $s_3=6 / 7$ 这一灰度级, 所以有 $656+329=985$ 个像素取这个值。同样, 有 $245+122+81=448$ 个像素取 $s_4=1$ 这个新灰度值。用 $n=4096$ 来除上述这些 $n_k$ 值, 便可得到新的直方图, 如图4-10(c)所示。

由上例可见,利用累积分布函数作为灰度变换函数,经变换后得到的新直方图虽然不很平坦,但毕竟比原始图像的直方图平坦的多,而且其动态范围也大大地扩展了。因此,这种方法对于对比度较弱的图像进行处理是很有效的。图4-11是经直方图均衡化后的Lena图像和其直方图。 因为直方图是近似的概率密度函数,所以用离散灰度级作变换一般得不到完全平坦的结果。另外,从上例可以看出,变换后的灰度级减少了,这种现象叫做“简并”现象。由于简并现象的存在,处理后的灰度级总是要减少的,这是像素灰度有限的必然结果。由于上述原因,数字图像的直方图均衡只是近似的

image-20230423192537615

图4-11 经直方图均衡化后的Lena图像及直方图 (a) 经直方图均衡化后的Lena图像; (b) 均衡化后的Lena图像的直方图

直方图均衡化

期望输出图像概率密度函数是均匀分布的 $$ \begin{aligned} & P_g(g)=\frac{1}{g_{\text {max }}-g_{\text {min }}}, g_{\text {min }} \leq g \leq g_{\text {max }}, \ & \int_{g_{\text {min }}}^g \frac{1}{g_{\text {max }}-g_{\text {min }}} d \omega=C(f), \ & \frac{1}{g_{\text {max }}-g_{\text {min }}}\left(g-g_{\text {min }}\right)=C(f), \ & g=\left[g_{\text {max }}-g_{\text {min }}\right] C(f)+g_{\text {min }} \end{aligned} $$

image-20230423192617430

直方图均衡化算法

  1. 列出原始图像灰度级 $f_j, j=0,1, \ldots, k, \ldots, L-1$;
  2. 统计各灰度级的象素数目 $, n_j, j=0,1,, \ldots, k, \ldots, L-1$;
  3. 计算原始图像直方图 $P_f\left(f_j\right)=n_j / n, n$ 为原始图像总的象素数目;
  4. 计算累积分布函数 $c(f) ; \quad c(f)=\sum_{j=0}^m P_g\left(g_i\right), i=0,1, \ldots, m, \ldots, L-1$ 5.应用转移函数, 计算映射后的灰度级,

$$ g_i=I N T\left[\left(g_{\max }-g_{\text {min }}\right) c(f)+g_{\min }+0.5\right] $$

(0.5是为了取整后达到四舍五入的效果)

  1. 统计映射后各灰度级的象素数目 $n_i, i=0,1, \ldots, p-1$;
  2. 计算输出图像直方图 $P_g\left(g_i\right)=n_j / n, i=0,1, \ldots, p-1$;
  3. $f_j$$g_i$ 的映射关系, 修改原始图像灰度级, 获得直方图近似均匀分布的输出图像。

典型直方图变换的转移函数

$$ \begin{array}{|l|l|l|} \hline 输出图像灰度概率密度数学模型 &数学模型表达式 &转移函数\ \hline 均匀(Uniform)& P_g(g)=\dfrac{1}{g_{\max }-g_{\min }}, g_{\max } \geq g \geq g_{\min} & g=\left[g_{\max }-g_{\min }\right] C(f)+g_{\min }\ \hline 指数(Exponitial)&P_g(g)=\alpha \exp \left[-\alpha\left(g-g_{\min }\right)\right], g \geq g_{\min }&g=g_{\min }-\dfrac{1}{\alpha} \ln [1-C(f)]\ \hline 瑞利(Raleigh)&P_g(g)=\dfrac{g-g_{\min }}{\alpha^2} \exp \left[-\dfrac{\left(g-g_{\min }\right)^2}{2 \alpha^2}\right], g \geq g_{\min }&g=g_{\min }+\left{2 \alpha^2 \ln \left[\dfrac{1}{1-C(f)}\right]\right}^{1 / 2}\ \hline 双曲线(立方根,Hyperbolic)&P_g(g)=\frac{1}{3}\left[\dfrac{g^{-2 / 3}}{g^{1 / 3}\max -{g^{1 / 3}\min }}\right], g_\max \geq g \geq g_{\min }&g=\left{\left[g^{1 / 3}{ }{\max }-g^{1 / 3}{ }{\min }\right] C(f)+g_{\min }^{1 / 3}\right}^3\ \hline 双曲线(对数,Hyperbolic)&P_g(g)=\dfrac{1}{g\left[\ln \left(g_{\max }\right)-\ln \left(g_{\min }\right)\right]}, g_{\max } \leq g \leq g_{\min }&g=g_{\min }\left[\dfrac{g_{\max }}{g_{\min }}\right]^{C(f)}\ \hline \end{array} $$

4.3 灰度变换

4.2.1 灰度线性变换

假定原图像 $f(x, y)$ 的灰度范围为 $[a, b]$, 希望变换后图像 $g(x$, $y)$ 的灰度范围扩展至 $[c, d]$, 则线性变换可表示为 $$ g(x, y)=\frac{d-c}{b-a}[f(x, y)-a]+c\tag{4-12} $$ 此式用图4-12表示。若图像灰度在 $0 \sim M_f$ 范围内, 其中大部分像素的灰度级分布在区间 $[a, b]$, 很小部分的灰度级超出了此区间, 为改善增强的效果, 可令: $$ g(x, y)= \begin{cases}c & 0 \leq f(x, y)<a \ \dfrac{d-c}{b-a}[f(x, y)-a]+c & a \leq f(x, y)<b \ d & b \leq f(x, y)<M_f\end{cases}\tag{4-13} $$ 有时为了保持$f(x, y)$灰度低端和高端值不变,可以采用如式(4-14)所示的形式 $$ g(x, y)=\left{\begin{array}{l} \dfrac{d-c}{b-a}[f(x, y)-a]+c \quad a \leq f(x, y) \leq b \ f(x, y) \end{array}\right.\tag{4-14} $$ 式中的$a$、$ b$、$ c$、$ d$这些分割点可根据用户的不同需要来确定。

例如, 当 $a=50, b=80, c=20, d-c=120$ 时, 即采用下式: $$ g(x, y)= \begin{cases}\dfrac{120}{30}[f(x, y)-50]+20 & 50 \leq f(x, y)<80 \ f(x, y) & \text { 其他 }\end{cases} $$ image-20230423194034171

图4-13 线性灰度变换

(a) 原始图像; (b) 灰度变换后的图像

线性灰度变换效果比较

image-20230423194050112

4.2.2 分段线性变换

为了突出感兴趣的目标或灰度区间,相对抑制那些不感兴趣的灰度区域,可采用分段线性变换。常用的三段线性变换法如图4-14所示,其数学表达式如下: $$ g(x, y)= \begin{cases}\dfrac{c}{a} f(x, y) & 0 \leq f(x, y)<a \ \dfrac{d-c}{d-a}[f(x, y)-a]+c & 0 \leq f(x, y)<b \ \dfrac{M_g-d}{M_f-b}[f(x, y)-b]+d & b \leq f(x, y)<M_f\end{cases}\tag{4-16} $$ image-20230423194202739

图4-15 常见的几种非线性变换函数

式(4-16)对灰度区间$ [0, a]$和$[b, M_f]$加以压缩,对灰度区间[a, b]进行扩展。通过细心调整折线拐点的位置及控制分段直线的斜率,可对任一灰度区间进行扩展或压缩。这种变换适用于在黑色或白色附近有噪声干扰的情况。例如, 照片中的划痕, 由于变换后在[0, a]以及[b, $M_f$ ]之间的灰度受到压缩, 因而使污斑得到减弱。

4.2.3 非线性变换

image-20230423194246479

S型函数

$$ \begin{aligned} & T(x)=\frac{L}{2}\left{1+\frac{1}{\sin \left(\alpha \frac{\pi}{2}\right)} \sin \left[\alpha \pi\left(\frac{r}{L}-\frac{1}{2}\right)\right]\right} & 0<\alpha<1 \\ & ~T(x)=\frac{D_m}{2}\left{1+\frac{1}{\tan \left(\alpha \frac{\pi}{2}\right)} \tan \left[\alpha \pi\left(\frac{x}{D_m}-\frac{1}{2}\right)\right]\right} & 0<\alpha<1 \end{aligned} $$

指数变换

$$ s=c y^\gamma $$

$c$$y$ 为常数。 当 $\gamma&lt;1$ 时, 具有扩展低亮度区、压缩高亮度区的效果; 当 $\gamma&gt;1$ 时, 性能正好相反,具有压缩低亮度区、扩展高亮度区的效果。 而当 $c=\gamma=1$ 时, 输出正比与输入, 图像无变化。 “ $\gamma$ 校正"(伽马校正),CRT显示器

显示器校正

image-20230423194412359

图6-23 不同$\gamma$的$s=cr^{\gamma}$曲线及图像变换结果

组合灰度变换效果与直方图变图像增强与平滑化

组合灰度变换效果与直方图变化

image-20230423194605299

分段线性变换

image-20230423194641215

4.4 图像噪声

重点:均值滤波法s

4.4.1 概述

噪声可以理解为“妨碍人们感觉器官对所接收的信源信息理解的因素” 。例如,一幅黑白图像,其亮度分布假定为$f (x, y)$, 那么对其起干扰作用的亮度分布$R(x, y)$便称为图像噪声。噪声在理论上可以定义为“不可预测, 只能用概率统计方法来认识的随机误差” 。因此,将图像噪声看成是多维随机过程是合适的, 描述噪声的方法完全可以借用随机过程及其概率分布函数和概率密度函数。但在很多情况下,这种描述方法很复杂,甚至不可能,而且实际应用往往也不必要,通常是用其数字特征,即均值方差、相关函数等进行处理。

目前,大多数数字图像系统中,输入光图像都是采用先冻结再扫描方式将多维图像变成一维电信号,再对其进行处理、存储、 传输等加工,最后往往还要再组成多维图像信号。图像噪声同样也受到这样的分解和合成,在这些过程中电气系统和外界影响将使得图像噪声的精确分析变得十分复杂。另一方面, 对图像信息的认识和理解是由人的视觉系统所决定的。不同的图像噪声,人的感觉(理解)程度是不同的,这就是所谓人的噪声视觉特性问题。该方面虽早已进行研究,但终因人的视觉系统本身未搞清楚而未获得解决。尽管如此,图像噪声在数字图像处理技术中的重要性却愈加明显。例如,高放大倍数遥感图片的判读,X射线图像系统中的噪声去除等都已成为不可缺少的技术。

4.4.2 图像噪声分类

图像噪声按其产生的原因可分为外部噪声和内部噪声。

外部噪声是指系统外部干扰从电磁波或经电源传进系统内部而引起的噪声,如电气设备、天体放电现象等引起的噪声。 主要外部干扰如下:

(1) 由光和电的基本性质所引起的噪声。

(2) 电器的机械运动产生的噪声。如, 各种接头因抖动引起的电流变化所产生的噪声;磁头、磁带抖动引起的抖动噪声等。

(3) 元器件材料本身引起的噪声。如, 磁带、磁盘表面缺陷所产生的噪声。

(4) 系统内部设备电路所引起的噪声。如, 电源系统引入的交流噪声,偏转系统和箝位电路引起的噪声等。

图像噪声从统计特性可分为平稳噪声和非平稳噪声两种。 统计特性不随时间变化的噪声称为平稳噪声;统计特性随时间变化的噪声称为非平稳噪声

另外, 按噪声和信号之间的关系可分为加性噪声和乘性噪声。假定信号为 $S(t)$, 噪声为 $n(t)$, 如果混合叠加波形是 $S(t)+n(t)$ 形式, 则称其为加性噪声; 如果叠加波形为 $S(t)[1+n(t)]$ 形式, 则称其为乘性噪声。为了分析处理方便, 往往将乘性噪声近似认为加性噪声, 而且总是假定信号和噪声是互相独立的。

4.4.3 图像系统噪声特点

如图4-16是一幅含有噪声的图像, 由此可知图像中的噪声有以下三个特点。

  1. 噪声在图像中的分布和大小不规则
  2. 噪声与图像之间具有相关性
  3. 噪声具有叠加性

image-20230423194824007

图4-16 有噪声的图像

4.5 去除噪声

改善降质图像的方法有两类:

一类是不考虑图像降质的原因, 只将图像中感兴趣的部分加以处理或突出有用的图像特征,故改善后的图像并不一定要去逼近原图像。这一类图像改善方法称为图像增强,主要目的是要提高图像的可懂度。

通过图像增强的方法提高图像质量,这一类方法称之为黑盒子办法.不去了解内部机理,只研究输入和输出

另一类方法是针对图像降质的具体原因,设法补偿降质因素,使改善后的图像尽可能地逼近原始图像。这类方法称为图像恢复或图像复原技术

图像增强处理的方法基本上可分为空间域法和频率法两大类。前者是在原图像上直接进行数据运算,对像素的灰度值进行处理。它又分为两类,一类是对图像作逐点运算,称为点运算;另一类是在与处理像点邻域有关的空间域上进行运算, 称为局部运算。频域法是在图像的变换域上进行处理, 增强感兴趣的频率分量, 然后进行反变换, 得到增强了的图像

图像恢复则被称之为白盒子办法,需要了解内部机理

4.5.1 模板操作和卷积运算

超重点:模板匹配法,掌握则保证85分

模板操作(前提:空域)是数字图像处理中常用的一种运算方式,图像的平滑、锐化以及后面将要讨论的细化、边缘检测等都要用到模板操作。例如, 有一种常见的平滑算法是将原图中的一个像素的灰度值和它周围邻近8个像素的灰度值相加,然后将求得的平均值作为新图像中该像素的灰度值。可用如下方法来表示该操作: $$ \frac{1}{9}\left[\begin{array}{ccc} 1 & 1 & 1 \ 1 & 1 * & 1 \ 1 & 1 & 1 \end{array}\right] $$

星号是代替的意思

上式有点类似于矩阵,通常称之为模板(Template),带星号的数据表示该元素为中心元素,即这个元素是将要处理的元素。如果模板为 $$ \frac{1}{9}\left[\begin{array}{ccc} 1 * & 1 & 1 \ 1 & 1 * & 1 \ 1 & 1 & 1 \end{array}\right] $$ 则该操作的含义是:将原图中一个像素的灰度值和它右下相邻近的8个像素值相加,然后将求得的平均值作为新图像中该像素的灰度值

模板操作实现了一种邻域运算,即某个像素点的结果不仅和本像素灰度有关,而且和其邻域点的值有关。模板运算的数学含义是卷积(或互相关)运算。 卷积是一种用途很广的算法,可用卷积来完成各种处理变换, 图4-17说明了卷积的处理过程

image-20230328165413026

图4-17 卷积运算示意图

卷积运算中的卷积核就是模板运算中的模板,卷积就是作加权求和的过程。邻域中的每个像素(假定邻域为3×3大小,卷积核大小与邻域相同),分别与卷积核中的每一个元素相乘,乘积求和所得结果即为中心像素的新值。卷积核中的元素称作加权系数(亦称为卷积系数),卷积核中的系数大小及排列顺序, 决定了对图像进行区处理的类型。改变卷积核中的加权系数, 会影响到总和的数值与符号, 从而影响到所求像素的新值。

在模板或卷积的加权运算中,还存在一些具体问题需要解决:首先是图像边界问题,当在图像上移动模板(卷积核)至图像的边界时,在原图像中找不到与卷积核中的加权系数相对应的9个像素,即卷积核悬挂在图像缓冲区的边界上, 这种现象在图像的上下左右四个边界上均会出现。例如,当模板为 $$ \frac{1}{9}\left[\begin{array}{ccc} 1 & 1 & 1 \ 1 & 1 * & 1 \ 1 & 1 & 1 \end{array}\right] $$ 设原图像为 $$ \left[\begin{array}{lllll} 1 & 1 & 1 & 1 & 1 \ 2 & 2 & 2 & 2 & 2 \ 3 & 3 & 3 & 3 & 3 \ 4 & 4 & 4 & 4 & 4 \end{array}\right] $$ 经过模板操作后的图像为 $$ \left[\begin{array}{ccccc}

  • & - & - & - & - \
  • & 2 & 2 & 2 & - \
  • & 3 & 3 & 3 & - \
  • & - & - & - & - \end{array}\right] $$ “-”表示无法进行模板操作的像素点。

解决这个问题可以采用两种简单方法:一种方法是忽略图像边界数据, 另一种方法是在图像四周复制原图像边界像素的值,从而使卷积核悬挂在图像四周时可以进行正常的计算。实际应用中,多采用第一种方法。

其次,是计算出来的像素值的动态范围问题, 对此可简单地将其值置为0或255即可。

4.5.2 邻域平均法

邻域平均法是一种利用Box模板对图像进行模板操作(卷积运算)的图像平滑方法,所谓Box模板是指模板中所有系数都取相同值的模板, 常用的3×3和5×5模板如下: $$ \frac{1}{9}\left[\begin{array}{ccc} 1 & 1 & 1 \ 1 & 1 * & 1 \ 1 & 1 & 1 \end{array}\right] $$

$$ \frac{1}{25}\left[\begin{array}{ccc} 1 & 1 & 1 & 1&1 \\ 1 & 1 & 1 & 1&1 \\ 1 & 1 & 1* & 1&1 \\ 1 & 1 & 1 & 1&1 \\ 1 & 1 & 1 & 1&1 \\ \end{array}\right] $$

Box模板对当前像素及其相邻的的像素点都一视同仁,统一进行平均处理, 这样就可以滤去图像中的噪声。例如,用3×3 Box模板对一幅数字图像处理结果,如图4-18所示(图中计算结果按四舍五入进行了调整,对边界像素不进行处理)。 $$ \begin{array}{|c|c|c|c|c|} \hline 1 & 2 & 1 & 4 & 3 \ \hline 1 & 2 & 2 & 3 & 4 \ \hline 5 & 7 & 6 & 8 & 9 \ \hline 5 & 7 & 6 & 8 & 8 \ \hline 5 & 6 & 7 & 8 & 9 \ \hline \end{array}\Rightarrow \begin{array}{|c|c|c|c|c|} \hline 1 & 2 & 1 & 4 & 3 \ \hline 1 & 3 & 4 & 4 & 4 \ \hline 5 & 4 & 5 & 6 & 9 \ \hline 5 & 6 & 7 & 8 & 8 \ \hline 5 & 6 & 7 & 8 & 9 \ \hline \end{array} $$ 图4-18 3×3Box模板平滑处理示意图

邻域平均法的数学含义可用下式表示: $$ g(x, y)=\frac{1}{M} \sum_{(i, j) \in s} f(i, j)\tag{4-17} $$ 式中: $x, y=0,1, \ldots, N-1 ; S$ 是以 $(x, y)$ 为中心的邻域的集合, $M$$S$ 内的点数。

邻域平均法的思想是通过一点和邻域内像素点求平均来去除突变的像素点,从而滤掉一定的噪声,其主要优点是算法简单,计算速度快, 但其代价是会造成图像一定程度上的模糊。

采用邻域平均法对图4-19(a)中的图像进行处理后的结果如图4-19(b)所示。可以看出经过邻域平均法处理后,虽然图像的噪声得到了抑制,但图像变得相对模糊了。

image-20230328165858329

图4-19 图像的领域平均法

(a) 原始图像; (b) 邻域平均后的结果

邻域平均法的平滑效果与所采用邻域的半径(模板大小)有关。半径愈大, 则图像的模糊程度越大, 因此, 减少图像的模糊是图像平滑处理研究的主要问题之一。

为解决邻域平均法造成图像模糊的问题,可采用阈值法、K邻点平均法、梯度倒数加权平滑法、 最大均匀性平滑法、小斜面模型平滑法等。它们讨论的重点都在于如何选择邻域的大小、 形状和方向,如何选择参加平均的点数以及邻域各点的权重系数等。有关这些方法请参阅相关参考文献。

具有边界保持的平滑滤波器

  1. 设计思想:

    前面的处理结果可知,经过平滑滤波处理之后,图像会变得模糊。分析原因,在图像上的景物之所以可以辨认清楚是因为目标物之间存在边界。

    这类滤波器是增加了一个判别当前像素点平滑时,选择非边界点进行。

  2. K近邻均值(中值)滤波器

    1. 以待处理像素为中心,作一个m*m的作用模板。

    2. 在模板中,选择K个与待处理像素的灰度差为最小的像素。

    3. 将这K个像素的灰度均值(中值)替换掉原来的像素值。

例:3*3模板,k=5

image-20230328170058303

4.5.3 中值滤波

中值滤波是一种非线性信号处理方法,与其对应的中值滤波器也就是一种非线性滤波器。中值滤波器于1971年提出并应用在一维信号时间序列分析中, 后来被二维图像信号处理技术所引用。它在一定条件下,可以克服线性滤波器(如邻域平滑滤波等)所带来的图像细节模糊,而且对滤除脉冲干扰及图像扫描噪声最为有效。在实际运算过程中并不需要图像的统计特性,这也带来不少方便。但是对一些细节多,特别是点、线、尖顶细节多的图像不宜采用中值滤波。

由于中值滤波是一种非线性运算,对随机输入信号的严格数学分析比较复杂,下面采用直观方法简要介绍中值滤波的原理。

中值滤波原理

中值滤波就是用一个奇数点的移动窗口, 将窗口中心点的值用窗口内各点的中值代替。假设窗口内有五点,其值为80、 90、 200、 110和120, 那么此窗口内各点的中值即为110。

设有一个一维序列 $f_1, f_2, \ldots, f_n$, 取窗口长度 (点数) 为 $m(m$ 为奇数), 对其进行中值滤波, 就是从输入序列中相继抽出 $m$ 个数 $f_{i-v}, \ldots, f_{i-1}, f_i, f_{i+1}, \ldots, f_{i+v}$ (其中 $f_i$ 为窗口中心点值, $v=(m-1) / 2$ ), 再将这 $m$ 个点按其数值大小排序, 取其序号为中心点的那个数作为滤波输出。用数学公式表示为 $$ y_i=\operatorname{Med}\left{f_{i-v}, \Lambda, f_i, \Lambda f_{i+v}\right} \quad i \in N, v=\frac{m-1}{2}\tag{4-18} $$ 例如, 有一序列 ${0,3,4,0,7}$, 重新排序后为 ${0,0,3,4,7}$, 则Med ${0,0,3,4,7}=3$ 。此列若用平滑滤波, 窗口也是取 5 , 那么平滑滤波输出为 $(0+3+4+0+7) / 5=2.8$ 。 图4-20使用内含 5 个像素的窗口对离散阶跃函数、斜坡函数、 脉冲函数以及三角形函数进行中值滤波和平均值滤波的示例。 左边一列为原波形, 中间为平均值滤波结果, 右边为中值滤波结果。可以看出, 中值滤波器不影响阶跃函数和斜坡函数。周期小于 $m / 2$ (窗口之半) 的脉冲受到抑制, 另外三角函数的顶部变平。

image-20230423195415339

图4-20 中值滤波和平均值滤波比较 (a) 阶跃;(b)斜坡;(c)单脉冲;(d)双脉冲;(e)三脉冲; (f)三角波

二维中值滤波可由下式表示: $$ y_{i j}={Med}A\left{f{i j}\right}\tag{4-19} $$ 式中: $A$ 为窗口; $\left{f_{i i}\right}$ 为二维数据序列。

二维中值滤波的窗口形状和尺寸对滤波效果影响较大,不同的图像内容和不同的应用要求,往往采用不同的窗口形状和尺寸。常用的二维中值滤波窗口有线状、方形、圆形、十字形以及圆环形等。窗口尺寸一般先用3×3,再取5×5逐渐增大, 直到滤波效果满意为止。就一般经验来讲,对于有缓变的较长轮廓线物体的图像,采用方形或圆形窗口为宜。对于包含有尖顶物体的图像, 用十字形窗口, 而窗口大小则以不超过图像中最小有效物体的尺寸为宜。如果图像中点、线、尖角细节较多,则不宜采用中值滤波。

中值滤波主要特性

对某些输入信号中值滤波的不变性

对某些特定的输入信号, 如在窗口内单调增加或单调减少的序列, 中值滤波输出信号仍保持输入信号不变, 即: $f_{i-}$ ${ }n \leq \ldots \leq f_i \leq \ldots \leq f{i+n}$ 或 $f_{i-n} \geq \ldots \geq f_i \geq \ldots \geq f_{i+n}$, 则 $\left{y_i\right}=\left{f_i\right}$ 一维中值滤波这种不变性可以从图4-20中 $(a)$$(b)$ 上看出来。 二维中值滤波的不变性如图4-21所示。它不但与输入信号有关, 而且还与窗口形状有关。图4-22给出了几种二维窗口及与之对应的最小尺寸的不变输入图形。一般与窗口对顶角连线垂直的边缘线保持不变性。利用这个特点, 可以使中值滤波既能去除图像中的噪声, 又能保持图像中一些物体的边缘。


image-20230328171136466

图4-21 中值滤波不变性示例

(a) 原始图像;(b) 中值滤波输出

image-20230328171200177

图4-22 中值滤波几种常用窗口及其相应的不变图形

对于一些周期性的数据序列,中值滤波也存在着不变性。

例如,下列一维周期性二值序列

$$ \left{f_{i}\right}=\ldots,+1,+1,-1,-1,+1,+1,-1,-1, \ldots $$

若设窗口长度为 9 , 则中值滤波对此序列保持不变性。对于二维周期序列不变性, 如周期网状结构图案, 分析起来就更复杂了, 可以通过试验改变窗口形状和尺寸来获取。

中值滤波去噪声性能

对于零均值正态分布的噪声输入, 中值滤波输出的噪声方差 $\sigma^{2}{ }_{\text {med }}$ 近似为

$$ \sigma_{m e d}^{2}=\frac{1}{4 m f^{2}(\bar{m})} \approx \frac{\sigma_{i}^{2}}{m+\frac{\pi}{2}-1} \cdot \frac{\pi}{2}\tag{4-20} $$

式中: $\sigma_{i}^{2}$ 为输入噪声功率 (方差) , $m$ 为中值滤波窗口长度 (点数) , $\bar{m}$ 为输入噪声均值, $f(\bar{m})$ 为输入噪声密度函数。

而均值滤波的输出噪声方差 $\sigma^{2}{ }{0}$ 为 $$ \sigma{0}^{2}=\frac{1}{m} \sigma_{i}^{2}\tag{4-21} $$

比较公式(4-20)和(4-21), 可以看出, 中值滤波的输出与输入噪声的密度分布有关。对随机噪声的抑制能力, 中值滤波比平均值滤波要差一些。但对脉冲干扰, 特别是脉冲宽度小于 $m$ / 2、相距较远的窄脉冲干扰,中值滤波的效果较好。

中值滤波的频谱特性

$G$ 为输入信号频谱, $F$ 为输出信号频谱, 定义中值滤波的频率响应特性为

$$ H=\left|\frac{G}{F}\right|\tag{4-22} $$

试验表明, $H$$G$ 的关系曲线如图4-23所示。由图可见, 中值滤波频谱特性起伏不大, 其均值比较平坦。可以认为信号经中值滤波后, 频谱基本不变。这一特点对设计和使用中值滤波器很有意义。

图4-23 $H$$G$ 的关系曲线

在计算机上制作了一组用于平滑实验的图像, 如图4-24所示。图(d)和(e)分别是采用 $3 \times 3$ 窗口算术平均平滑去除噪声后的图像,图(f)和(g)分别是采用 $5 \times 5$ 十字中值滤波后的图像。显然, 算术平均平滑对含有高斯噪声的图像有效; 而中值滤波对含有椒盐噪声图像的去噪声效果较好。

image-20230328171345472

图4-24 噪声平滑实验图像 (a)Lena原图; (b) 高斯噪声; (c) 椒盐噪声; (d) 对(b)平均平滑; (e) 对(c)平均平滑; (f) 对(b) $5 \times 5$ 中值滤波; (g) 对(c) $5 \times 5$ 中值滤波

复合型中值滤波

对一些内容复杂的图像, 可以使用复合型中值滤波。如, 中值滤波线性组合、高阶中值滤波组合、加权中值滤波以及迭代中值滤波等。

\section{(1)中值滤波的线性组合是将几种窗口尺寸大小和形状不同} 的中值滤波器复合使用, 只要各窗口都与中心对称, 滤波输出可保持几个方向上的边缘跳变, 而且跳变幅度可调节。其线性组合方程如下:

$$ Y_{i j}=\max {k}\left[\mathop{Med}\limits{A_k}\left{f_{i j}\right}\right]\tag{4-23} $$

式中: $A_{k}$ 为窗口。

(2) 高阶中值滤波组合如下式所示:

$$ Y_{i j}=\sum_{k=1}^{N} a_{k} \mathop{Med}\limits_{A_k}\left{f_{i j}\right}\tag{4-24} $$

式中: $a_{k}$ 为不同中值滤波的系数, $A_{k}$ 为窗口。

这种中值滤波可以使输入图像中任意方向的细线条保持不变。例如,可选择图4-25中的4种线状窗口 $A_{1} \sim A_{4}$, 用式(4-24)组合的中值滤波, 可以使输入图像中各种方向的线条保持不变, 而且又有一定的噪声平滑性能。

image-20230328171523765

图4-25 几种线性窗口

(3)其他类型的中值滤波: 为了在一定的条件下尽可能去除噪声, 又有效保持图像细节, 可以对中值滤波器参数进行修正, 如加权中值滤波, 也就是对输入窗口进行加权。也可以是对中值滤波器的使用方法进行变化, 保证滤波的效果, 还可以和其他滤波器联合使用。

4.5.4 其他去噪技术

空间域低通滤波

从信号频谱角度来看, 信号的缓慢变化部分在频率域属于

低频部分, 而信号的迅速变化部分在频率域是高频部分。对图像来说, 它的边缘以及噪声干扰的频率分量都处于频率域较高的部分, 因此, 可以采用低通滤波的方法来去除噪声。而频域的滤波又很容易从空间域的卷积来实现, 为此只要适当设计空间域的单位冲激响应矩阵, 就可以达到滤除噪声的效果。下面是几种用于噪声平滑低通卷积模板。

$$ \frac{1}{4}\left[\begin{array}{ccc}0 & 1 & 0 \ 1 & 0 * & 1 \ 0 & 1 & 0\end{array}\right]\quad\quad \frac{1}{8}\left[\begin{array}{ccc}0 & 1 & 0 \ 1 & 0 * & 1 \ 1 & 1 & 1\end{array}\right]\quad\quad \frac{1}{10}\left[\begin{array}{ccc}1 & 1 & 1 \ 1 & 2 * & 1 \ 1 & 1 & 1\end{array}\right] $$ 这些模板中引入了加权系数,以区分邻域中不同位置像素对输出像素值的影响, 常称其为加权模板。与邻域平均法中采用的 Box模板相比较可知,Box模板并没有考虑邻域中各点位置的影响,对于所有的邻点都一视同仁,所以其平滑的效果并不理想。

另外, 二维Gaussian离散模板也是一种常用的低通卷积模板。 由于Gaussian函数有着一些良好的特性, 对二维连续Gaussian分布经采样、量化, 并使模板归一化, 便可得到二维Gaussian离散模板。3×3二维Gaussian模板如下: $$ \frac{1}{16}\left[\begin{array}{ccc} 1 & 2 & 1 \ 2 & 4^{*} & 2 \ 1 & 2 & 1 \end{array}\right] $$

可见Gaussian离散模板也是一种加权模板, 并且它是按二维正态分布进行加权的。

s频率域低通滤波

利用卷积定理, 可以把式(4-10)写成以下形式:

$$ G(u, v)=H(u, v) F(u, v)\tag{4-25} $$

式中: $F(u, v)$ 是含噪声图像的傅立叶变换, $G(u, v)$ 是平滑后图像的傅立叶变换, $H(u, v)$ 是低通滤波器传递函数。

利用 $H(u, v)$ 使 $F(u, v)$ 的高频分量得到衰减, 得到 $G(u, v)$ 后再经过反变换就得到所希望的图像 $g(x, y)$ 。低通滤波平滑图像的系统框图如图4-26所示。

image-20230423195944002

图4-26 低通滤波平滑图像的系统框图

​ 常用的频率域低通滤波器有:理想低通滤波器(ILPF)、 巴特沃思低通滤波器(BLPF)、指数低通滤波器(ELPF)。这三种频率域低通滤波器的频率特性比较如图4-27所示。

image-20230328171755708

图4-27 ILPF、BLPF、ELPF特征曲线 (a) ILPF特征曲线; (b) BLPF特征曲线; (c) ELPF特征曲线

多幅图像平均法

一幅有噪声的图像 $f(x, y)$, 可以看作是由原始无噪声图像 $g(x, y)$} 和噪声 $n(x, y)$ 叠加而成(加性噪声),即

$$ f(x, y)=g(x, y)+n(x, y)\tag{4-26} $$

若叠加在图像上的噪声 $n(x, y)$ 是非相关、具有零均值的随机噪声时,那么,把针对一目标物(景物)在相同条件下, 把作M次重复摄取的图像相加, 取平均值作为输出图像,便可对图像中的噪声进行平滑。

多幅图像的平均输出图像的期望值是无噪声的理想图像 $g(x$,y)。需要注意的是对多幅图像平均, 要求多幅图像之间相互对准,而大多数图像要做到严格对准是相当困难的。多幅图像取平均处理常用于摄像机的视频图像中, 用以减少电视摄像机光电摄像管或CCD器件所引起的噪声。这是对同一景物连续摄取多幅图像并数字化, 再对多幅图像平均, 一般选用8幅图像取平均。这种方法在实际应用中的难点是如何把多幅图像配准起来, 以便使相应的像素能正确地对应排列。

要明确的一点是:没有一种十全十美的去除噪声方法, 应根据处理目的,通过实际试验选择合适的方法。

4.6 图像锐化

需要使用高通滤波器,此时变化大的点为感兴趣的点

所谓锐化算法,就是利用图像的一阶微分信息,得到想要的信息

4.6.1 微分法

梯度法

对于图像函数 $f(x, y)$, 它在点 $(x, y)$ 处的梯度是一个矢量, 定义为 $$ G[f(x, y)]=\left[\begin{array}{l} \frac{\partial f}{\partial x} \ \frac{\partial f}{\partial y} \end{array}\right] \tag{4-27} $$ 梯度的两个重要性质是:

(1) 梯度的方向在函数$f(x, y)$最大变化率的方向上。

(2) 梯度的幅度用 $G[f(x, y)]$ 表示, 并由下式算出: $$ G[f(x, y)]=\left[\left(\frac{\partial f}{\partial x}\right)^2+\left(\frac{\partial f}{\partial y}\right)^2\right]^{1 / 2}\tag{4-28} $$ 由上式可知, 梯度的数值就是 $f(x, y)$ 在其最大变化率方向上的单位距离所增加的量。对于数字图像而言, 式(4-28)可以近似为 $$ G[f(x, y)]=\left{[f(i, j)-f(i+1, j)]^2+[f(i, j)-f(i, j+1)]^2\right}^{1 / 2}\tag{4-29} $$ 式中: 各像素的位置如图4-28(a)所示。

image-20230317102022338

图4-28 求梯度的两种差分运算

式(4-29)可简化成为 $$ G[f(x, y)]=\left|f(i, j)-f(i+1, j)\right|+\left|f(i, j)-f(i, j+1)\right| \tag{4-30} $$ 以上梯度法又称为水平垂直差分法。另一种梯度法叫做罗伯特梯度法 (Robert Gradient), 它是一种交叉差分计算法, 如图4-28(b)所示。其数学表达式为 $$ G[f(x, y)]=\left{[f(i, j)-f(i+1, j+1)]^2+[f(i+1, j)-f(i, j+1)]^2\right}^{1 / 2}\tag{4-31} $$ 同样可近似为 $$ G[f(x, y)]=\mid\left[f(i, j)-f\left(i+1, j^{+1)}|+| f(i+1, j)-f(i, j+1)\right.\right.\tag{4-32} $$ 以上两种梯度近似算法在图像的最后一行和最后一列的各像素的梯度无法求得,一般就用前一行和前一列的梯度值近似代替。 由梯度的计算可知,图像中灰度变化较大的边缘区域其梯度值大,在灰度变化平缓的区域其梯度值较小,而在灰度均匀区域其梯度值为零。图4-29(b)是采用水平垂直差分法对图4-29(a)锐化的结果,锐化后仅留下灰度值急剧变化的边沿处的点

image-20230317102159818

图4-29 图像梯度锐化结果

(a) 二值图像; (b) 梯度运算结果

当梯度计算完之后,可以根据需要生成不同的梯度增强图像。 第一种是使各点的灰度$g(x, y)$等于该点的梯度幅度,即 $$ g(x, y)=G[f(x, y)]\tag{4-33} $$ 此法的缺点是增强的图像仅显示灰度变化比较陡的边缘轮廓, 而灰度变化平缓的区域则呈黑色

第二种增强的图像是使 $$ g(x, y)= \begin{cases}G[f(x, y)] & G[f(x, y)] \geq T \ f(x, y) & \text {其他 }\end{cases}\tag{4-34} $$ 式中: $T$ 是一个非负的阈值, 适当选取 $T$, 即可使明显的边缘轮廓得到突出, 又不会破坏原灰度变化比较平缓的背景。

第三种增强图像是使 $$ g(x, y)= \begin{cases}L_G & G[f(x, y)] \geq T \ f(x, y) & \text {其他 }\end{cases}\tag{4-35} $$ 式中: $T$ 是根据需要指定的一个灰度级, 它将明显边缘用一固定的灰度级 $L_G$ 来实现。

第四种增强图像是使 $$ g(x, y)= \begin{cases}G[f(x, y)] & G[f(x, y)] \geq T \ L_G & \text { 其他 }\end{cases}\tag{4-36} $$ 此法将背景用一个固定灰度级 $L_G$ 来实现, 便于研究边缘灰度的变化。

第五种增强图像是使 $$ g(x, y)= \begin{cases}L_G & G[f(x, y)] \geq T \ L_B & \text { 其他 }\end{cases}\tag{4-37} $$ 此法将背景和边缘用二值图像表示, 便于研究边缘所在位置。

Sobel算子

重点,很常用,目前国际上占据主流

采用梯度微分锐化图像,同时会使噪声、条纹等得到增强, Sobel算子则在一定程度上克服了这个问题。Sobel算子法的基本原理是: 按式(4-38)计算3×3窗口(如图4-30所示)的灰度, 将其作为变换后图像f(i, j)的灰度。

image-20230317102637380

图4-30 Sobel算子图像坐标 $$ g=\sqrt{S_x^2+S_y^2}\tag{4-38} $$ 式中: $$ \begin{aligned} S_x & =[f(i+1, j-1)+2 f(i+1, j)+f(i+1, j+1)] \ & -[f(i-1, j-1)+2 f(i-1, j)+f(i-1, j+1)] \ S_y & =[f(i-1, j+1)+2 f(i, j+1)+f(i+1, j+1)] \ & -[f(i-1, j-1)+2 f(i, j-1)+f(i+1, j-1)] \end{aligned} $$ 为简化计算, 可用 $\left.g=\left|S_{x}\right|\right]+\left|S_{y}\right|$ 来代替式(4-38), 从而得到锐化后的图像。Sobel算子不像普通梯度算子那样用两个像素的差值, 这就导致了以下两个优点: (1) 由于引入了平均因素, 因而对图像中的随机噪声有一定的平滑作用。 (2) 由于它是相隔两行或两列之差分, 故边缘两侧元素得到了增强, 边缘显得粗而亮。

Sobel锐化算法 $$ d_x=\left[\begin{array}{rrr} -1 & 0 & 1 \ -2 & 0 & 2 \ -1 & 0 & 1 \end{array}\right]&d_y=\left[\begin{array}{ccc} -1 & -2 & -1 \ 0 & 0 & 0 \ 1 & 2 & 1 \end{array}\right]\tag{4-39} $$

$$

$$

Roberts梯度锐化效果图例

image-20230317102827205

为简化计算, 可用 $g=\left|S_{x}\right|+\left|S_{y}\right|$ 来代替式(4-38), 从而得到锐化后的图像。Sobel算子不像普通梯度算子那样用两个像素的差值, 这就导致了以下两个优点: (1) 由于引入了平均因素, 因而对图像中的随机噪声有一定的平滑作用。 (2) 由于它是相隔两行或两列之差分, 故边缘两侧元素得到了增强, 边缘显得粗而亮。

Roberts梯度锐化效果图例

image-20230423200429260

Priwitt锐化算法 $$ g(i, j)=\left{d_x^2(i, j)+d_y^2(i, j)\right}^{\frac{1}{2}}\tag{4-40}\ $$ $$ d_x=\left[\begin{array}{ccc} -1 & 0 & 1 \ -1 & 0 & 1 \ -1 & 0 & 1 \end{array}\right] &d_y=\left[\begin{array}{ccc} -1 & -1 & -1 \ 0 & 0 & 0 \ 1 & 1 & 1 \end{array}\right] $$

表4-2 常用的梯度算子 $$ \begin{array}{|l|l|l|ll|} \hline \text { 算子名称 } & H_1 & H_2 & \text { 特 } & \text { 点 } \ \hline \text { Roberts } & {\left[\begin{array}{cc} 0 * & 1 \ -1 & 0 \end{array}\right]} & {\left[\begin{array}{cc} 1 * & 0 \ 0 & -1 \end{array}\right]} & \text { 边缘定位准,但对噪声敏感 } \ \hline \text { Prewitt } & {\left[\begin{array}{ccc} -1 & 0 & 1 \ -1 & 0 * & 1 \ -1 & 0 & 1 \end{array}\right]} & {\left[\begin{array}{ccc} -1 & -1 & 1 \ 0 & 0 * & 0 \ 1 & 1 & 1 \end{array}\right]} & \text { 平均、微分对噪声有抑制作用 } \ \hline \text {Sobel } & {\left[\begin{array}{ccc} -1 & 0 & 1 \ -2 & 0 * & 2 \ -1 & 0 & 1 \end{array}\right]} & {\left[\begin{array}{ccc} -1 & -1 & 1 \ 0 & 0 * & 0 \ 1 & 2 & 1 \end{array}\right]} & 加权平均边宽 \geqslant 2 像素 \ \hline \text { Krisch } & {\left[\begin{array}{ccc} 5 & 5 & 5 \ -3 & 0 * & -3 \ -3 & -3 & -3 \end{array}\right]} & {\left[\begin{array}{ccc} -3 & -3 & 5 \ -3 & 0 * & 5 \ -3 & -3 & 5 \end{array}\right]} & 对噪声有较好地抑制作用。该算子需求出 f(x, y) 在 8 个方向上的平均差分的最大值, 这里只给出了 2 个方向的模板 \ \hline \text { Isotropic Sobel } & {\left[\begin{array}{ccc} -1 & 0 & 1 \ -\sqrt{2} & 0 * & \sqrt{2} \ -1 & 0 & 1 \end{array}\right]} & {\left[\begin{array}{ccc} -1 & -1 & 1 \ 0 & 0 * & 0 \ 1 & \sqrt{2} & 1 \end{array}\right]} & \text { 权值反比于领点与中心点的距离, 检测沿不同方向边缘时梯度幅度一致 } \ \hline \end{array} $$ Prewitt锐化效果图例

image-20230423201738704

其他

一、单方向的一阶梯度算法(浮雕效果)

水平方向的锐化

image-20230317103237686

一阶水平方向锐化效果

image-20230317103258348

  1. 垂直方向的锐化

image-20230317103357018

单方向一阶锐化效果图例

image-20230317103422861

image-20230317103451436

4.6.2 拉普拉斯运算

拉普拉斯算子是常用的边缘增强算子,拉普拉斯运算也是偏导数运算的线性组合运算,而且是一种各向同性(旋转不变性)的线性运算。拉普拉斯算子为 $$ \nabla^2 f=\frac{\partial^2 f}{\partial x^2}+\frac{\partial^2 f}{\partial y^2}\tag{4-41} $$ 如果图像的模糊是由扩散现象引起的 (如胶片颗粒化学扩散等), 则锐化后的图像 $g$ 为 $$ g=f-k \nabla^2 f\tag{4-42} $$ 式中: $f 、 g$ 分别为锐化前后的图像, $k$ 为与扩散效应有关的系数。

式(4-42)表示模糊图像经拉普拉斯算子锐化以后得到新图像 $g \circ k$ 的选择要合理, 太大会使图像中的轮廓边缘产生过冲; $k$ 太小,锐化不明显。 对数字图像来讲, $f(x, y)$ 的二阶偏导数可表示为 $$ \begin{aligned} \frac{\partial^2 f(x, y)}{\partial x^2} & =\nabla_x f(i+1, j)-\nabla_x f(i, j) \ & =[f(i+1, j)-f(i, j)]-[f(i, j)-f(i-1), j] \ & =f(i+1, j)+f(i-1, j)-2 f(i, j) \ \end{aligned} $$ $$ \frac{\partial^2 f(x, y)}{\partial y^2} =f(i, j+1)+f(i, j-1)-2 f(i, j)\tag{4-43} $$

为此, 拉普拉斯算子 $\nabla^2 f$ 为 $$ \begin{aligned} \nabla^2 f & =\frac{\partial^2 f(x, y)}{\partial x^2}+\frac{\partial^2 f(x, y)}{\partial y^2} \ & =f(i+1, j)+f(i-1, j)+f(i, j+1)+f(i, j-1)-4 f(i, j) \ & =-5\left{f(i, j)-\frac{1}{5}[f(i+1, j)+f(i-1, j)+f(i, j+1)\right. \ & +f(i, j-1)+f(i, j)]} \end{aligned}\tag{4-44} $$ 可见, 数字图像在 $(i, j)$ 点的拉普拉斯算子, 可以由 $(i, j)$ 点灰度值减去该点邻域平均灰度值来求得。当 $k=1$ 时, 拉普拉斯锐化后的图像为 $$ \begin{aligned} g(i, j) & =f(i, j)-\nabla^2 f(i, j) \ & =5 f(i, j)-f(i+1, j)-f(i-1, j)-f(i, j+1)-f(i, j-1) \end{aligned}\tag{4-45} $$

Laplacian锐化算子

写成模板系数形式形式即为Laplacian算子: $$ H_1=\left[\begin{array}{ccc} 0 & -1 & 0 \ -1 & 4 & -1 \ 0 & -1 & 0 \end{array}\right] $$

Laplacian变形算子

$$ H_2=\left[\begin{array}{ccc} -1 & -1 & -1 \ -1 & 8 & -1 \ -1 & -1 & -1 \end{array}\right] \quad H_3=\left[\begin{array}{ccc} 1 & -2 & 1 \ -2 & 4 & -2 \ 1 & -2 & 1 \end{array}\right] \quad H_4=\left[\begin{array}{ccc} 0 & -1 & 0 \ -1 & 5 & -1 \ 0 & -1 & 0 \end{array}\right] $$ 例:

设有 $1 \times n$ 的数字图像 $f(i, j)$, 其各点的灰度如下:

$$ \ldots, 0,0,0,1,2,3,4,5,5,5,5,5,5,6,6,6,6,6,6,3,3,3,3,3 $$

计算 $\nabla^2 f$ 及锐化后的各点灰度值 $g$ (设 $\left.k=1\right)$

  • 首先按式 (4-44)计算各点的 $\nabla^2 f$ 。例如:

第3点:

$$ \quad \nabla^2 f=-3\left[0-\frac{1}{3}(0+0+1)\right]=1 $$

第8点:

$$ \quad \nabla^2 f=-3\left[5-\frac{1}{3}(4+4+5)\right]=-1 $$

各点拉普拉斯算子如下:

$$ \ldots, 0,0,1,0,0,0,0,-1,0,0,0,0,1,-1,0,0,0,0,-3,3,0,0,0, \ldots $$

  • 按式(4-42)计算 $g=f-2 f$ 。例如:

$$ \begin{aligned} &第3点: g=f-\nabla^2 f=0-1=-1\\ &\vdots \\ &第8点: g=5-(-1)=6 \end{aligned} $$

锐化后各点的灰度值如下:

$$ \ldots, 0,0,-1,1,2,3,4,6,5,5,5,5,4,7,6,6,6,6,9,0,3,3,3, \ldots $$

从上例可以看出,在灰度级斜坡底部(如第3点)和界线的低灰度级侧(如第13、20点)形成下冲。在灰度级斜坡顶部(如第8点)和界线的高灰度级侧(如第14、 19点)形成上冲。在灰度级平坦区域(如第912点, 第1518点), 运算前后没有变化。如图4-31所示。

image-20230423202505266

图4-31 拉普拉斯锐化前、 后图像的灰度 (a) 原图像灰度; (b) 拉普拉斯锐化后图像的灰度

image-20230525145157358

图4-33 拉普拉斯锐化结果 (a) 二值图像; (b) 拉普拉斯运算结果

Laplacian锐化算子效果

image-20230317103909298

image-20230317104014112

Wallis算子 $$ g(i, j)=\log [f(i, j)]-\frac{1}{4} s \ s=[\log f(i-1, j)+\log f(i+1, j)+\log f(i, j-1)+\log f(i, j+1) $$ 它可以看作是校正了视觉的指数特性后的Laplacian运算。 $$ H_1=\left[\begin{array}{ccc} 0 & -1 & 0 \ -1 & 4 & -1 \ 0 & -1 & 0 \end{array}\right] \Rightarrow H_1=\left[\begin{array}{ccc} 0 & -\dfrac{1}{4} & 0 \ -\dfrac{1}{4} & 1 & -\dfrac{1}{4} \ 0 & -\dfrac{1}{4} & 0 \end{array}\right] $$ Wallis算法效果图

image-20230525145608343

image-20230317104139383

图4-34 3种频域高通滤波特性曲线 (a) IHPF特征曲线; (b) BHPF特征曲线; (c) EHPF特征曲线

4.7 图像的伪彩色处理

4.7.1 密度分割

密度分割是伪彩色处理技术中最简单的一种。设一幅灰度图像 $f(x, y)$, 在某一个灰度级 (如 $\left.f(x, y)=L_1\right)$ 上设置一个平行于 $x y$ 平面的切割平面, 如图4-35所示。灰度图像被切割成只有两个灰度级, 对切割平面以下的 (灰度级小于 $L_1$ ) 像素分配给一种颜色 (如蓝色), 对切割平面以上的像素分配给另一种颜色 (如红色)。这样切割结果就可以将灰度图像变为只有两个颜色的伪彩色图像。

若将灰度图像级用 $M$ 个切割平面去切割。就会得到 $M+1$ 个不同灰度级的区域 $S_1, S_2, \ldots, S_{M}, S_{M+1}$ 。对这 $M+1$ 个区域中的像素人为分配给 $M+1$ 种不同颜色, 就可以得到具有 $M+1$ 种颜色的伪彩色图像, 如图4-36所示。密度分割伪彩色处理的优点是简单易行, 便于用软件或硬件实现。还可以扩大它的用途, 如计算图像中某灰度级面积等。

image-20230317110351336

图4-35 密度分割示意图

image-20230317110407463

图4-36 多灰度伪彩色分割示意图

4.7.2 灰度级彩色变换

这种伪彩色处理技术(在遥感技术中常称为假彩色合成方法),可以将灰度图像变为具有多种颜色渐变的连续彩色图像, 实际图像的连续伪彩色变换如图4-37所示。其变换过程为:将灰度图像送入具有不同变换特性的红、绿、蓝3个变换器,再将3 个变换器的不同输出分别送到彩色显像管的红、绿、蓝电子枪。 同一灰度由3个变换器对其实施不同变换, 而使3个变换器输出不同,从而在彩色显像管里合成某种色彩。可见,不同大小灰度级一定可以合成不同色彩。

image-20230317110436524

图4-37 伪彩色变换

从图中可见, 若 $f(x, y)=0$, 则 $I_B(x, y)=L, I_R(x, y)=I_{G}(x, y)=0$, 从而显示蓝色。同样, 若 $f(x, y)=L / 2$, 则 $I_G(x, y)=L, I_R(x, y)=I_G(x$, $y)=0$, 从而显示绿色。若 $f(x, y)=L$, 则 $I_R(x, y)=L, I_B(x, y)=I_G(x$, $y)=0$, 从而显示红色。 因此不难理解, 若灰度图像 $f(x, y)$ 灰度级在 $0 \sim L$ 之间变化, $I_R 、 I_B 、 I_G$ 会有不同输出, 从而合成不同的彩色图像。

4.7.3 滤波法

这是一种在频率域进行伪彩色处理的技术,与上面不同的是输出图像的伪彩色与图像的灰度级无关, 而是取决于图像中不同空间频率域成分。如, 为了突出图像中高频成分(图像细节)将其变为红色,只要将红色通道滤波器设计成高通特性即可。而且可以结合其他处理方法,在附加处理中实施(如直方图修正等),使其彩色对比度更强。如果要抑制图像中某种频率成分,可以设计一个带阻滤波器。其过程如图4-38所示,从3个不同频率的滤波器输出的信号再经过傅立叶反变换,可以对其做进一步的处理,如直方图均衡化。最后把它们作为三基色分别加到彩色显像管的红、绿、蓝显示通道,从而实现频率域的伪彩色处理。

image-20230425192923253

图4-38 频率域伪彩色增强处理

4.8常见滤波器

4.8.1频域低通滤波

ILPF

$$ \begin{aligned} & H(u, v)=\left{\begin{array}{l} 1, \text { ifD }(u, v) \leq D_0 \\ 0, \text { ifD }(u, v)>D_0 \end{array}\right. \\ & D(u, v)=\sqrt{u^2+v^2} \end{aligned} $$

image-20230425193206897

图4-39理想低通滤波器(a)变换函数透视图;(b)图像方式显示的滤波器(c)滤波器的剖面图

截止频率与所保留能量的关系 $$ E_T=\sum_{u=0}^{N-1} \sum_{v=0}^{N-1} E(u, v) \ E(u, v)=R^2(u, v)+I^2(u, v) \ \beta=100\left[\sum_u \sum_v E(u, v) / E_T\right] \ $$ image-20230425193301672

图4-40低通滤波器的能量与$D_0$的关系(a)尺寸为$500\times500$像素的图像(b)图像的Fourier谱。叠加的圆环分别具有5、15、30、80、230像素的半径,分别包含的能量为90.0%、94.6%、96.4%、98.0%、99.55%

IHPF的图像处理效果

image-20230425193359163

图4-41 图像4-39(a)经过理想低通滤波器的结果, 对应的半径分别为15、30、80,滤去的能量为5.4%、3.6%、2%。

振铃效应

$$ G(u,v)=H(u,v)F(u,v)\\ g(x,y)=h(x,y)*f(x,y) $$

image-20230425193518474

图4-42(a)半径为5的脉冲图像 (b)相应的空间滤波器

(c)空域的5个脉冲(d)滤波结果

滤波器1

image-20230425193601972

滤波器2

image-20230425193616333

巴特沃思滤波器Butterworth

BLPF: $$ H(u, v)=\frac{1}{1+\left[D(u, v) / D_0\right]^{2 n}} $$ image-20230425193703918

图4-43 BLPF特性曲线

BLPF的滤波结果

image-20230425193721411

图4-44 二阶BLPF滤波的结果(a)原图像, (b)半径15,(b)半径30,(d)半径80

指数滤波器或高斯滤波器

ELPF: $$ H(u, v)=e^{-\left[D(u, v) / D_0\right]^n} $$ image-20230425193751790

图4-45高斯滤波器(a)传递函数的透视图;(b)图像形式;(c)径向剖面图

ELPF的滤波效果

image-20230425193807471

图4-46 高斯滤波器的处理效果, 从左到右依次为原图、半径15、30、80

4.8.2高通滤波器

image-20230425193827184

图4-47 三种典型的高通滤波器,从上到下依次为理想高通滤波器、巴特沃思高通滤波器、指数滤波器。

高通滤波器定义

IHPF: $$ H(u, v)=\left{\begin{array}{l} 0, \text { if } D(u, v) \leq D_{0 ;} \ 1, \text { if } D(u, v)>D_0 . \end{array} \quad D(u, v)=\sqrt{u^2+v^2}\right.\ $$ BHPF: $$ H(u, v)=\frac{1}{1+\left[D_0 / D(u, v)\right]^{2 n}} $$ EHPF: $$ H(u, v)=e^{-\left[D_0 / D(u, v)\right]^{n}} $$ THPF: $$ H(u, v)=\left{\begin{array}{lr} 0 & D(u, v)<D_0 \ \dfrac{\left[D(u, v)-D_1\right]}{D_0-D_1} & D_0 \leq D(u, v) \leq D_1 \ 1 & D(U, v)>D_1 \end{array}\right. $$ HPF比较

image-20230425194000220

高通滤波器效果

image-20230425194028009

图4-48 理想高通滤波效果,从左到右依次为$D_0$=15,30,80。越小,振铃效应越明显

image-20230425194049850

图4-49 BHPF效果,从左到右依次为$D_0$=15,30,80。比IHPF的结果平滑得多

image-20230425194146423

图5-50高斯高通滤波效果,从左到右依次为$D_0$=15,30,80。

高频加强

$$ \begin{aligned} & H^{\prime}(u, v)=H(u, v)+c, 0<c<1 . \\ & G^{\prime}(u, v)=F(u, v) H^{\prime}(u, v)=F(u, v) H(u, v)+c F(u, v) \end{aligned} $$

image-20230425194231541

图4-51 X光片原图;高通滤波效果;高频增强效果;直方图均衡化效果

4.8.3带通(带阻)滤波器

抑制以点 $(u 0, v 0)$ 为中心, $D 0$ 为半径的邻域中所有频率的理想带阻滤波器 (IBPF) 的转移函数为: $$ H(u, v)=\left{\begin{array}{lll} 0, & \text { if } & D(u, v) \leq D_0 \ 1, & \text { if } & D(u, v)>D_0 \end{array}\right.\ D(u, v)=\left{\left(u-u_0\right)^2+\left(v-v_0\right)^2\right}^{1 / 2} $$ image-20230425194334938

不围绕原点的情况

$$ H(u,v)=\begin{cases} 0 &&if,,D_1(u,v)\leq D_0 ,,or,, D_2(u,v)\leq D_0\\ 1&& else \end{cases} $$

$$ \begin{aligned} &\begin{aligned} & D_1(u, v)=\left{\left(u-u_0\right)^2+\left(v-v_0\right)^2\right}^{1 / 2} \\ & \end{aligned}\\ &D_2(u, v)=\left{\left(u+u_0\right)^2+\left(v+v_0\right)^2\right}^{1 / 2} \end{aligned} $$

围绕原点的模型

为消去围绕原点的一个频带,必须考虑周期性和对称性。一个径向对称(放射对称)的理想带阻滤波器 $$ H(u, v)=\left{\begin{array}{rr} 1, & \text { if } D(u, v)<D_0-\frac{W}{2} \ 0, & \text { if } D_0-\frac{W}{2} \leq D(u, v) \leq D_0+\frac{W}{2} \ 1, & \text { if } D(u, v)>D_0+\frac{W}{2} \end{array}\right. $$ image-20230425194440894

巴特沃思带阻滤波器用于去噪声

image-20230425194459895

(a)被正弦噪声污染的图像 (b)图像(a)的频谱 (c)巴特沃思带阻滤波器 (d)图像(a)的滤波效果

频域滤波小结

物理可实现的低通滤波器的传递函数必须是连续的。如果处处可导,滤波后的图像就不会出现振铃效应;

频域低通滤波器的效果是去除图像中的高频噪声,去噪声的能力与滤波器的形式以及截止频率有关;

频域低通滤波器在滤除噪声的同时会造成图像模糊,图像模糊的程度与截止频率有关;

高通滤波器可以增强图像的边缘,是实现图像边缘检测的基础;

高通滤波器基本上抑制了图像中的平滑信息,因此如要在保持图像基本信息的基础上实现锐化,就需要采用高频加强技术;

带通和带阻滤波器在消除特定类型的噪声上有很好的效果,它们也是彩色图像增强的技术基础;

同态滤波是一类利用图像的照明反射模型,同时实现亮度动态范围压缩和图像对比度增强的技术

4.9 编程实例

4.7.1 编程绘制灰度直方图

1. 创建灰度直方图数据

//******************************************* 
//函数名称: CreateHistogram()
//参数说明: nX1、 nY1、 nX2、 nY2为需创建直方图图像区域坐标, 默认为- // 	          1。pData、pPalette、 nWidthBytes为位图数据指针、调色板指
// 	        针、 图像行字节宽度
//返回值: int型数组指针*pBuffer。在数组中依次存储亮度、红、绿、蓝分量
//	   的直方图数据, 每个分量各占256个元素[FK)]
//基本功能:创建指定的图像对象的直方图数据。如果未给定坐标,则针对整	       个图像
//**************************************** 
 int *CPointPro:: CreateHistogram(int nX1, int nY1, int nX2, int nY2,
         unsigned char *pData, RGBQUAD *pPalette, int nWidthBytes)
{
  //图像指针为空, 无法操作返回
  if(m_pDibObject == NULL) return(FALSE); 
  //分配直方图数据缓存区(数组)
  int *pBuffer = new int256 * 4]; 
  //分配直方图数据缓存区失败
  if(pBuffer == NULL) return( NULL ); 
  //直方图数据缓存区清零 
memset(pBuffer, 0, (256 * 4) * sizeof(int)); 
  //变量定义
  DWORD dwGray;
  int x, y;
  unsigned char *pTemp, ucRed, ucGreen, ucBlue;
  //图像的高度
  int nHeight = m_pDibObject->GetHeight();
  switch(m_pDibObject->GetNumBits())
  {
  case 24: 
    for(y = nY1; y < nY2; y++)
    {
      	//数据指针定位到图像数据起始位置
	pTemp = pData;
      	//数据指针定位到图像数据每行的起始零位置 
      pTemp += (( nHeight - 1 - y) * nWidthBytes);
      //数据指针定位到图像数据每行的起始nX1位置
      pTemp += (nX1 * 3);
      for(x=nX1; x<=nX2; x++)
      {
        //获取像素颜色的三原色
        ucRed  = pTemp[x * 3 + 2];
        ucGreen = pTemp[x * 3 + 1];
        ucBlue  = pTemp[x * 3];
        //按关系L=0.3R+0.59G+0.11B, 得到亮度值
	dwGray  = (DWORD) (ucRed * 30 + ucGreen * 59 + ucBlue * 11) / 100;
        dwGray &= 0x000000ff; 
       //亮度直方图数据
        pBuffer[dwGray]++;
        //红色直方图数据
        pBuffer[256 + ucRed]++; 
        //绿色直方图数据
        pBuffer[512 + ucGreen]++;
        //蓝色直方图数据
        pBuffer[768 + ucBlue]++;
        //数据指针加3
        pTemp += 3;
               }
        }
        break;
     }
     return( pBuffer );
}  

2.定制灰度直方图对话框

定制一个直方图对话框并创建一个CDlgIntensity类, 在它的构造函数中初始化成员变量(例如, 初始化直方图数据指针m_pnHistogram = NULL), 添加WM_PAINT的消息映射函数OnPaint(), 然后在消息映射函数中绘制直方图数据。

3. 显示直方图

void CDipView:: OnViewHist() 
{
  CDipDoc* pDoc = GetDocument();
  ASSERT_VALID(pDoc);
  //判断当前是否有图像对象
  if( pDoc->m_pDibObject = = NULL ) return;
  //在点处理CPointPro类中创建用来绘制直方图的数据
  CPointPro PointOperation( pDoc->m_pDibObject );
  int *pHistogram = PointOperation.GetHistogram();
  //生成一个对话框CDlgIntensity类的实例
  CDlgIntensity HistDlg; 
    //将直方图数据传给CDlgIntensity的成员变量m_pnHistogram
  if( pHistogram != NULL )
  {
    //设置直方图数据指针
    HistDlg.m_pnHistogram = pHistogram;
    //设置当前像素值为0的像素数
    HistDlg.m_nCurrentPiexsNum = pHistogram[0];
    //设置是否为256级灰度图像
    HistDlg.m_bIsGray256 = PointOperation.IsGray256();
  }
  //显示对话框
  if ( HistDlg.DoModal() != IDOK) return;
  delete [] pHistogram;
}

4.7.2 中值滤波的VC++编程实现

下面以3×3矩形窗口为例来说明8位灰度图像的中值滤波的主要算法。

for(y = nY1; y < nY2; y++)
{
  // pOldBits和pNewBits分别为指向源图像和新图像的图像数据指针
  pOldTemp = pOldBits; 
  pOldTemp += ((nHeight-1-y-1)*nWidthBytes); 
  pOldTemp += nX1-1; 
  pNewTemp = pNewBits; 
  pNewTemp += ((nHeight-1-y-1)*nWidthBytes); 
  pNewTemp += nX1; 
  k = 0; 
  for(x = nX1; x<nX2; x++)
  { 
    //取得3×3邻域像素的灰度值, 存放在dwMedianList数组中
    for(i = 0; i < 3; i++)
    {
      Data = pOldTemp[i+k];
      dwMedianList[i] = (pOldPalette[Data].rgbRed+pOldPalette[Data].rgbGreen +
        pOldPalette[Data].rgbBlue)/3;
      //数组nLoc存放3×3邻域像素在位图数据中的位置
      nLoc[i]=(i+k);
    } 
for( i=0; i<3; i++ )
{
      Data = pOldTemp[i+k+nWidthBytes];
      dwMedianList[i+3] = (pOldPalette[Data].rgbRed +
        pOldPalette[Data].rgbGreen +pOldPalette[Data].rgbBlue) / 3;
      nLoc[i+3] = ( i + k ) + nWidthBytes;
}
for( i=0; i<3; i++ )
{
      Data = pOldTemp[i+k+nWidthBytes*2];
      dwMedianList[i+6] = (pOldPalette[Data].rgbRed +
            pOldPalette[Data].rgbGreen +pOldPalette[Data].rgbBlue ) / 3;
      nLoc[i+6] = ( i + k ) + nWidthBytes * 2;
}
    //对3×3邻域像素的灰度值进行排序
    for( i=1; i < 9; i++ )
    {
      for( j=0; j<i; j++ )
      {
        if( dwMedianList[i] < dwMedianList[j])
        {
          DWORD dwTmp;
          int nTmp;
          dwTmp = dwMedianList[i];
          nTmp = nLoc[i];
          for( k = i; k>j; k- -)
          {
		dwMedianList[k] = dwMedianList[k-1
           		 nLoc[k] = nLoc[k-1];
          	  }
                dwMedianList[j] = dwTmp;
                nLoc[j]= nTmp;
                break;
            }
        }
    } 
    // 取得3×3邻域排序结果的中间值在原图像中的索引, 并给新图赋值
    Data = pOldTemp[nLoc[4]]; 
    *pNewTemp++ = (unsigned char) m_pDibObject->GetNearestIndex
         (pNewPalette[Data].rgbRed,  pNewPalette[Data].rgbGreen, 
         pNewPalette[Data].rgbBlue, pNewPalette, nNumColors );
    k++;
}
}  

4.7.3 图像平滑的VC++编程实现

  1. 设计模板对话框 图像平滑实际上可以看作是卷积运算, 因此图像平滑的关键在于模板窗口的设计。本例首先定制了如图4-39所示的对话框, 并创建了与之对应的CDlgSmooth类。除了可以使用平均模板和高斯模板外, 还可自定义模板, 模板的大小被限制在5×5内。

image-20230425194920303

图4-39 图像平滑模板设置对话框

   利用类向导为对话框的相应控件添加成员变量,如模板类型、模板高度、模板宽度、模板中心元素X坐标、模板中心元素Y坐标对应的成员变量分别定义为整型变量m_intType、m_iTempH、m_iTempW、m_iTempMX、m_iTempMY, 而模板系数对应浮点型的变量m_fTempC,另外还定义了一个模板元素数组指针m_fpArray。
    为3个单选按钮分别定义一个单击事件的消息映射,在映射函数中给模板元素m_fpArray及其他成员变量赋值。 

    2. 获取图像数据
    对每个像素的邻域用选定的模板作卷积运算, 就可以得到平滑后的图像。下面给出了卷积运算的源代码。
    卷积函数中的各参数的含义如下: fpData为被卷积的数组, fpKernel为卷积核数组,fCoef为卷积系数,nSize为卷积尺寸, nResult为卷积结果。 
void CAreaPro:: Convolution(float *fpData, float *fpKernel, float fCoef,
				     int nSize, unsigned char *nResult)
{
  int i;
  float Sum = 0, fResult;
  //计算卷积
  for(i = 0; i < nSize; i++)
  {
    Sum += fpData[i] * fpKernel[i];
  }
  //执行卷积后的结果
  fResult = Sum / fCoef; 
  //求绝对值
  fResult = (float)fabs(fResult);
  //判断是否超过255
  if(fResult > 255.0 )
  {
    // 直接赋值为255
    fResult = 255.0;
          }
  //对卷积结果四舍五入, 并转换成unsigned char类型作为最后返回结果
  *nResult = (unsigned char) (fResult + 0.5);
} 

​ 接着就是对原图像的每个像素作卷积运算,运算后的结果赋予新图像,其中图像边界像素不作处理,代码实现在CAreaPro:: Smoothing()中。