第五章 图像恢复和重建 CHAPTER 5 IMAGE RESTORATION and RECONSTRUCTION 第五章 图像恢复和重建 CHAPTER 5 IMAGE RESTORATION and RECONSTRUCTION §1 退化的数学模型和对角化 §2 无约束恢复 §3 有约束恢复 §4 几何失真校正 §5 投影重建 版权所有， 1997 (c) Dale Carnegie & Associates, Inc.

§5.1 退化的数学模型和对角化 退化：图像质量的降低； 失真可看作是退化，校正是恢复； 投影可看作是退化（三维到二维平面），重建是恢复； §5.1.1 简单的通用图像退化模型 　 n（x,y） f（x,y） H + g（x,y） 模型化：一个作用在f（x,y）上的系统H与一个加性噪声n（x,y） 的联合作用，导致产生退化图像g（x,y）。 　　假设已知n（x,y）的统计特性（或先求出），图像复原就是已知g（x,y）求f（x,y）的问题 （近似过程）。 g（x,y）= H [f（x,y）] + n（x,y） ； 已知 退化 解 噪声

§5.1.2 退化模型的计算 一、1D情况： 设f（x）和h（x）均匀采样后放入尺寸为A、B的数组； §5.1.2 退化模型的计算 一、1D情况： 设f（x）和h（x）均匀采样后放入尺寸为A、B的数组； f（x），x=0,1,…,A-1； h（x），x=0,1,…,B-1； 为了避免卷积的周期重叠，取M  A+B-1； 将f（x）和h（x）用零扩展补齐；(M周期)  ge（x）=  fe（m）h（x-m） ； x, m=0,1,…,M-1； 用矩阵形式表示的卷积为 g（x）= H f ，展开得 ge（0） = he（0） he（-1） … he（-M+1） fe（0） ge（1） = he（1） he（0） … he（-M+2） fe（1） ge（M-1）= he（M-1） he（M-2）… he（0） fe（M-1）

§5.1.2 退化模型的计算（续1） 式中， he（0） he（-1） … he（1） H矩阵是一个循环矩阵，其特点有：行列首尾相接；环形 §5.1.2 退化模型的计算（续1） 式中， he（0） he（-1） … he（1） H = he（1） he（0） … he（2） he（M-1） he（M-2）… he（0） H矩阵是一个循环矩阵，其特点有：行列首尾相接；环形 循环矩阵相加后仍是循环矩阵，相乘后仍是循环矩阵；

§5.1.2 退化模型的计算（续2） 二、2D情况（1D的直接推广，0x  M-1, 0  y  N-1 ） §5.1.2 退化模型的计算（续2） 二、2D情况（1D的直接推广，0x  M-1, 0  y  N-1 ） he(x,y)= h(x,y) ，0x  A-1, 0  y  B-1 ，否则he(x,y)= 0； fe(x,y)= f (x,y) ，0x  C-1, 0  y  D-1 ，否则fe(x,y)= 0；  ge (x,y) =   fe (m,n) he (x-m,y-n) ； m=0,1,…,M-1； n=0,1,…,N-1； 用矩阵形式表示的卷积为 g = H f ，展开得 ge（0） = H0 HM-1 … H1 fe（0） ge（1） = H1 H0 … H2 fe（1） ge（M-1）= HM-1 HM-2 … H0 fe（MN-1） 其中每个块Hi 是由扩展函数he (x,y)的第i行而来，是块循环矩阵。

§5.1.3 循环矩阵对角化 一、1D对角化 设循环矩阵H的本征值和特征向量分别为： §5.1.3 循环矩阵对角化 一、1D对角化 设循环矩阵H的本征值和特征向量分别为： (k) = [ 1 exp[j2/M *k] …… exp[j2/M* (M-1)k]T； (k) = he(0)+he(M-1)exp[j2/M•k]+……+ he(1)exp[j2/M • (M-1) k] ； 将H的M个特征向量组成一个M*M的矩阵W： W=[(0) (1) (2) …… (M-1)] H = WDW-1； 式中D为一个对角矩阵，其元素正是H本征值， 即 D (k,k)= (k)。

§5.1.3 循环矩阵对角化(续1) 设对4*4的循环矩阵C对角化。 §5.1.3 循环矩阵对角化(续1) 设对4*4的循环矩阵C对角化。 设(k)为整数1的M次根，k为根的序号，即 (k)M = 1 = lne； 因(k)=exp[j2k / M]；例中M=4；即(k)4 = 1； (k)=exp[j2k/4]； 规定标量(k)= c0+ (k) c1 + (k)2 c2 +(k)3 c3 = ∑ ci exp[j2k/4*i]； 根据线性代数分析，k=0,1,2,3时可得到四个本征值和相应的特征向量； (k) (k)=C (k)；组成(k)矩阵 W=[(0) (1) (2) (3)] 各元素 (k) =[ 1 exp (j2/4•k ) exp (j2/4•2k ) exp (j2/4•3k ) ]； 当循环矩阵已知本征值和特征向量后，可以进行对角化； 设对角化矩阵为D，则D= W-1CW；（ C等于H ）； D为对角阵，D (k,k)= (k) =MH (k) = W-1HW；

§5.1.3 循环矩阵对角化(续2) 二、2D对角化（块循环矩阵的对角化） 借助循环矩阵的结果进行推广； §5.1.3 循环矩阵对角化(续2) 二、2D对角化（块循环矩阵的对角化） 借助循环矩阵的结果进行推广； 设矩阵W的尺寸为MN*MN，每个元素为W(i,m)=exp(j2/M*im)WN ； WN为一个N*N矩阵，每个元素为WN(k,n)=exp(j2/N*kn) ； H = WDW-1 ；且HT = WD*W-1 ； D*为D的复共轭； ∴ D = W-1HW；可转化为对角阵。

§5.1.4 退化模型对角化的效果 一、1D情况： H=WDW-1 ；等式代入g=H f，并两边同左乘 W-1，得W-1 g= DW-1 f， W-1 g= DW-1 f 都是M维向量，令它们分别为G(k)和F(k)，得G(k)=DF(k)； 其中， G(k)= (1/M) ∑ ge (i) exp(j2/M*ki)是ge (x)的傅里叶变换； F(k)= (1/M) ∑ fe (i) exp(j2/M*ki)是fe (x)的傅里叶变换； 已知D (k,k)= (k) = ∑ he (i) exp(j2/M*ki)=MH (k) ， H (k)是he (i)的傅里叶变换， 所以G(k)= M× H (k) F(k)，式中 H (k) F(k)为fe (x) 与ge (x)的卷积， 可采用FFT计算。

§5.1.4 退化模型对角化的效果(续1) 二、2D情况： G(u,v)= (1/MN) ∑∑ ge (x,y) exp[-j2 (ux/M+vy/N)] ； F(u,v)= (1/MN) ∑∑ fe (x,y) exp[-j2 (ux/M+vy/N)] ； N(u,v)= (1/MN) ∑∑ ne (x,y) exp[-j2 (ux/M+vy/N)] ； H(u,v)= (1/MN) ∑∑ he (x,y) exp[-j2 (ux/M+vy/N)] ； D的MN个对角元素可用下式表示， D(k,i)= MN *H ([k/N]，k mod N) ；当i=k；否则为0； 其中[k/N]表示取最大整数，mod表示取余数； 所以G(u,v)= H(u,v) F(u,v) + N(u,v)；u,v=0,…,M-1

§5.1.5 有、无约束恢复 一、无约束恢复公式 f ‘ = (nTn) -1 HT g = H -1 (HT) –1 HT g = H-1 g ； 二、有约束恢复公式 若选取f ‘ 的一个线性操作符Q（变换矩阵）， 使得||Q f ‘ ||最小，采用拉格朗日乘数法解决，最小化下式的f ‘ ， L（ f ‘ ）= ||Qf ‘ ||2+ l（||g-H f ‘ || 2 -||n|| 2 ）， l 为拉格朗日乘数，该式为准则函数； 最小化得有约束恢复公式： f ‘ = [HT H + s QT Q ]-1 HT g；s = 1 / l ；

§5.2 无约束恢复 §5.2.1 逆滤波 一、无噪声时处理方法 §5.2 无约束恢复 §5.2.1 逆滤波 一、无噪声时处理方法 F ‘ (u,v) = G(u,v) / H(u,v) ；u,v=0,1,…,M-1 上式给定的恢复方法称为逆滤波，条件H –1 (u,v) 存在。 H(u,v)取零或很小时，计算上较困难。 二、有噪声时处理方法（G = H F + N，两端同除H） F ‘ (u,v) = F(u,v) + N(u,v) / H(u,v) ； H(u,v)取零或很小时，计算结果差距很大；一般情况下，逆滤波器并不正好是1 / H(u,v) ，记为 M(u,v)，称为恢复转移函数； M(u,v) = 1 / H(u,v) ；如u2+v2 ≤02 ， 0为去掉零点的值； M(u,v) = 1 ； 如u2+v2 > 02 ，

§5.2.1 逆滤波(续1) 图像退化和恢复模型 n(x,y) §5.2.1 逆滤波(续1) 图像退化和恢复模型 n(x,y) f (x,y) H(u,v) + g (x,y) M(u,v) f ‘ (x,y) 退化 恢复 改进：M (u,v) = k ； 如H(u,v) ≤ d ， d和k为小于1的常数； M(u,v) = 1/ H(u,v) ；其它， d越小越好 H (u,v) 已知，则可恢复原图像； 退化系统的转移函数可以用退化图像的傅里叶变换来近似，即 H (u,v) ≈ G (u,v) FFT 去零点 IFFT g (x,y) G(u,v) ≈ H(u,v) M(u,v) f ‘ (x,y) 点冲激下 F(u,v)=1

§5.2.2 运动模糊图像的恢复 运动模糊图像的恢复主要是消除匀速直线运动模糊，采用照相机摄影。 §5.2.2 运动模糊图像的恢复 运动模糊图像的恢复主要是消除匀速直线运动模糊，采用照相机摄影。 假设对平面匀速运动的景物采集一幅图像 f (x,y)，并设x0 (t) 和y0 (t) 分别是景物在x,y方向的运动分量； g(x,y) =∫0T f [x- x0 (t),y- y0 (t)] dt ；T是采集时间长度； g(x,y)即是由于运动而造成的模糊图像； 傅里叶变换 G(u,v)= F(u,v) H(u,v) ； 归结为无约束恢复类 H(u,v) = ∫0T exp[-j2(u x0 (t) +v y0 (t))] dt； 如果知道了运动分量x0 (t) 和y0 (t) ，就可解析得到传递函数H(u,v) 。