CN102096913A - 压缩感知框架下的多策略图像融合方法 - Google Patents

Abstract

本发明公开了一种压缩感知框架下的多策略图像融合方法，主要解决现有图像融合方法计算量大，时间复杂度高、存储空间大的问题，其实现过程是：输入原图像A和B，并将其分为大小为C×C的局部图像X1和X2；分别对X1和X2进行傅里叶变换，得到系数矩阵y1和y2；采用傅里叶系数低频全采的变密度观测模型，分别对y1和y2进行观测，得到观测向量f1和f2；根据f1和f2计算谐波系数H₁和H₂，以及频谱匹配度S；选取阈值T，计算加权系数；通过对加权系数、阈值和频谱匹配度的比较，计算融合后的观测向量f；对观测向量f进行Split Bregman重构算法迭代20次，最后得到所需融合图像。本发明与现有的融合方法相比具有计算复杂度低，融合效果好的优点，可用于视频跟踪、目标识别和计算机视觉。

Description

压缩感知框架下的多策略图像融合方法

技术领域

本发明属于图像处理技术领域，特别涉及图像的融合，可用于视频跟踪、目标识别和计算机视觉。

背景技术

图像融合是一种信息处理技术，即对多幅图像进行处理，以获得改善后的新图像。图像融合技术是研究如何加工、协同地利用多幅图像，并使不同图像信息相互补充，以获得对同一事物或目标地更客观、更本质认识的图像信息综合处理技术。由于可见光成像系统的聚焦范围有限和传感器的差异，因此在同一个场景中，聚焦良好的物体可呈现清晰的像，该物体前后一定距离外的所有目标都将呈现不同程度的模糊，而要得到一幅所有目标都清晰的像；由于传感器的差异造成的图像分辨率的差异，要想得到理想的视觉效果这都需要对图像进行融合。如何把从各种不同的传感器得到的图像融合起来，以便更充分地利用这些图像成为图像处理领域重要的研究课题之一。

传统的融合方法主要有：像素级融合、特征级融合和决策级融合。像素级图像融合是在严格配准的条件下，对多源图像直接进行图像融合处理，缺点是需要的数据量大，处理速度慢，实时性差。特征级图像融合是通过检测图像中的边缘、纹理和方向等特征信息提取细节特征再进行融合操作，融合效果比像素级融合有明显提高，缺点是实时性差。决策级融合是从图像中挖掘语义信息，再进行融合操作，优点是具有良好的实时性和容错性，缺点是预处理代价高，信息损失最多。

上述融合方法都是通过对图像的所有系数进行处理来进行图像融合，虽然具有到一定的融合效果，但还存在以下的不足：

(1)融合图像的数据量大，这导致图像存贮占用空间较大，不利于图像的压缩和传输；

(2)图像融合中数据的计算复杂度高，这使得图像融合过程时间较长。

发明内容

本发明的目的在于克服上述已有技术的不足，提出一种基于压缩感知的多策略图像融合方法，以减少图像传输数据量，降低图像融合的计算复杂度。

实现本发明目的的技术方案是：充分利用了图像的压缩感知理论，结合局部傅立叶域谐波的融合策略，通过少量的观测向量来进行图像融合，具体步骤如下：

(1)分别输入原图像A和原图像B，并将原图像A和原图像B分为大小为C×C的局部图像X1和X2，C×C取8×8或16×16；

(2)对局部图像X1进行傅里叶变换得到傅里叶系数矩阵y1，对局部图像X2进行傅里叶变换得到傅里叶系数矩阵y2；

(3)采用定密度和变密度相结合的观测模型，对傅里叶系数矩阵y1进行观测得到观测向量f1，对傅里叶系数矩阵y2进行观测得到观测向量f2；

(4)根据观测向量f1和观测向量f2计算两幅图像的谐波系数H₁、H₂以及频谱匹配度S：

$H_1=\sqrt{\frac{\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}{\left|{f1}_{i,j}\right|}^2}{\underset{m=1}{\overset{C}{\mathrm{\Σ}}}\underset{n=1}{\overset{C}{\mathrm{\Σ}}}\left|X1(m,n)\right|}}$

$H_2=\sqrt{\frac{\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}{\left|{f2}_{i,j}\right|}^2}{\underset{m=1}{\overset{C}{\mathrm{\Σ}}}\underset{n=1}{\overset{C}{\mathrm{\Σ}}}\left|X2(m,n)\right|}}$

$S=\frac{2\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}\left|{f1}_{i,j}\right|\×\left|{f2}_{i,j}\right|}{\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}{\left|{f1}_{i,j}\right|}^2+\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}{\left|{f2}_{i,j}\right|}^2}$

其中：H₁为图像A的谐波系数，X1(m，n)为局部图像X1在横坐标为m、纵坐标为n处的像素值，m∈[1，C]，n∈[1，C]，C为图像分块大小，f1_i，j为观测向量f1在横坐标位置为i、纵坐标位置为j处的值，i∈[1，C]，j∈[1，C]，H₂为图像B的谐波系数，X2(m，n)为局部图像X2在横坐标为m、纵坐标为n处的像素值，f2_i，j为观测向量f2在横坐标位置为i、纵坐标位置为j处的值；

(5)选取一个阈值T∈[0.6，1]，计算最大加权系数W_max和最小加权系数W_min：

W_max＝0.5-0.5×(1-S)/(1-T)

W_min＝1-W_max

其中S为步骤(4)中得到的频谱匹配度；

(6)根据步骤(4)中所得到的谐波系数H₁、H₂和频谱匹配度S，以及步骤(5)中所得到的最大加权系数W_max和最小加权系数W_min，计算融合后的观测向量f：

若S≥T且H1≥H2，则融合后的观测向量f＝W_max×f1+W_min×f2；

若S≥T且H1＜H2，则融合后的观测向量f＝W_min×f1+W_max×f2；

若S＜T且H1≥H2，则融合后的观测向量f＝f1；

若S＜T且H1＜H2，则融合后的观测向量f＝f2；

其中T为阈值，T∈[0.6，1]，f1为步骤(3)中得到的图像A的观测向量，f2为步骤(3)中得到的图像B的观测向量；

(7)根据步骤(6)中所得到的融合后的观测向量f，进行基于Wavelet-Curvelet框架的Split Bregman重构算法迭代20次，最后得到所需要的融合图像。

本发明与现有方法相比具有如下优点：

1.本发明由于在图像进行传输过程中，采用了傅里叶系数低频全采的变密度采样的模型，因而只传输少量的观测向量，而不是对图像数据全部进行传输，相对于传统的图像融合方法，减少了需要传输的数据量；

2.本发明由于在图像进行融合过程中，采用了傅里叶系数低频全采的变密度采样的模型，因而只需要对少量的观测向量进行融合，相对于传统的图像融合方法，由于融合数据量的减少，降低了融合的计算复杂度，缩短了融合时间；

3.本发明由于结合了多聚焦图像融合的特点，采取计算傅立叶域谐波系数和频谱匹配度来进行融合，因此对多聚焦图像有很好的融合效果。

附图说明

图1是本发明的流程图；

图2是本发明与现有平均值法对图像Clock1、Clock2的融合效果对比图；

图3是本发明与现有平均值法对图像Lab1、Lab2的融合效果对比图；

图4是本发明与现有平均值法对图像Pepsi1、Pepsi2的融合效果对比图；

图5是本发明与现有平均值法对图像Pepsi1、Pepsi2聚焦区域融合效果对比图。

具体实施方式

参照图1，本发明的具体实施过程如下：

步骤1，输入原图像A和原图像B，把原图像A和原图像B分为大小为C×C的局部图像X1和X2，C×C取8×8或16×16，本实例取16×16。

步骤2，对局部图像X1进行傅里叶变换得到傅里叶系数矩阵y1，对局部图像X2进行傅里叶变换到傅里叶系数矩阵y2。

步骤3，采用傅里叶系数低频全采的变密度观测模型，对傅里叶系数矩阵y1进行观测得到观测向量f1。

(3a)设定采样模型是取值仅为0或1的矩阵，将值为1的点作为采样点，根据输入图像A的大小设定矩阵B：若输入图像A的大小为m×m，则设矩阵B的大小为m×m且取值全为0，将坐标值为(0.5×m，0.5×m)的点a作为圆心，以0.3×0.3×m为半径的圆，将圆内所有点都作为采样点，将矩阵B中这些位置上的点值设为1；

(3b)按照

的采样密度对步骤(3a)圆外的采样点随机选取，取到的点值设为1，未取到的点值设为0，其中r为圆外任取一点b到圆心a的距离长度，R是以圆心a为起点通过点b，以图像边缘为终点的线段长度；

(3c)重复步骤(3a)和步骤(3b)z次，取z次的平均值作为最终的观测模型M，其中z取值为500到1000；

(3d)用步骤(3c)中得到的观测模型M点乘步骤2中得到的傅里叶系数矩阵y1，得到观测向量f1。

步骤4，采用傅里叶系数低频全采的变密度观测模型，对傅里叶系数矩阵y2进行观测，得到观测向量f2。

(4a)设定采样模型是取值仅为0或1的矩阵，将值为1的点作为采样点，根据输入图像A的大小设定矩阵B：若输入图像A的大小为m×m，则设矩阵B的大小为m×m且取值全为0，将坐标值为(0.5×m，0.5×m)的点a作为圆心，以0.3×0.3×m为半径的圆，将圆内所有点都作为采样点，将矩阵B中这些位置上的点值设为1；

(4b)按照的采样密度对步骤(4a)圆外的采样点随机选取，取到的点值设为1，未取到的点值设为0，其中r为圆外任取一点b到圆心a的距离长度，R是以圆心a为起点通过点b，以图像边缘为终点的线段长度；

(4c)重复步骤(4a)和步骤(4b)z次，取z次的平均值作为最终的观测模型M，其中z取值为500到1000；

(4d)用步骤(4c)中得到的观测模型M点乘步骤2中得到的对傅里叶系数矩阵y2进行观测，得到观测向量f2。

步骤5，根据观测向量f1和观测向量f2分别计算原图像A的谐波系数H₁和原图像B的谐波系数H₂，以及频谱匹配度S：

$H_1=\sqrt{\frac{\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}{\left|{f1}_{i,j}\right|}^2}{\underset{m=1}{\overset{C}{\mathrm{\Σ}}}\underset{n=1}{\overset{C}{\mathrm{\Σ}}}\left|X1(m,n)\right|}}$

$H_2=\sqrt{\frac{\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}{\left|{f2}_{i,j}\right|}^2}{\underset{m=1}{\overset{C}{\mathrm{\Σ}}}\underset{n=1}{\overset{C}{\mathrm{\Σ}}}\left|X2(m,n)\right|}}$

$S=\frac{2\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}\left|{f1}_{i,j}\right|\×\left|{f2}_{i,j}\right|}{\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}{\left|{f1}_{i,j}\right|}^2+\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}{\left|{f2}_{i,j}\right|}^2}$

其中：H₁为图像A的谐波系数，X1(m，n)为局部图像X1在横坐标为m、纵坐标为n处的像素值，m∈[1，C]，n∈[1，C]，C为图像分块大小，f1_i，j为观测向量f1在横坐标位置为i、纵坐标位置为j处的值，i∈[1，C]，j∈[1，C]，H₂为图像B的谐波系数，X2(m，n)为局部图像X2在横坐标为m、纵坐标为n处的像素值，f2_i，j为观测向量f2在横坐标位置为i、纵坐标位置为j处的值。

步骤6，选取一个阈值T∈[0.6，1]，计算最大加权系数W_max和最小加权系数W_min：

W_max＝0.5-0.5×(1-S)/(1-T)

W_min＝1-W_max

其中S为步骤5中得到的频谱匹配度。

步骤7，根据步骤5所得到的原图像A的谐波系数H₁、原图像B的谐波系数H₂，频谱匹配度S、最大加权系数W_max和最小加权系数W_min，计算融合后的观测向量f。

(7a)若S≥T且H1≥H2，则融合后观测向量f＝W_max×f1+W_min×f2；

(7b)若S≥T且H1＜H2，则融合后观测向量f＝W_min×f1+W_max×f2；

(7c)若S＜T且H1≥H2，则融合后观测向量f＝f1；

(7d)若S＜T且H1＜H2，则融合后观测向量f＝f2；

其中T为阈值，T∈[0.6，1]，f1为步骤3中得到的图像A的观测向量，f2为步骤4中得到的图像B的观测向量。

步骤8，根据步骤7所得到融合后的观测向量f，进行基于Wavelet-Curvelet框架的Split Bregman重构算法迭代20次，最后得到所需要的融合图像。

本发明的优点由以下仿真的数据和图像进一步说明。

1.仿真条件

1)选取六幅标准测试多聚焦待融合图像：Clock1、Clock2、Lab1、Lab2、Pepsi1、Pepsi2来测试融合效果，将本发明融合结果与现有的压缩感知框架下的平均值方法融合结果进行对比；

2)实验中将原图像Clock1和Clock2分为大小为16×16的块，阈值取T＝0.86；将原图像Pepsi1和Pepsi2分成大小为16×16的块，阈值取T＝0.9；将原图像Lab1和Lab2分成40×32的块，阈值取T＝0.95；

2.仿真内容和结果分析

1)本发明方法与现有平均值方法对图像融合的仿真实验

本实验的主要目的是对比本发明与现有平均值方法的图像融合效果，其实验数据采用六幅标准测试自然图像：Clock1、Clock2、Pepsi1、Pepsi2它们均是512×512维的，而Lab1和Lab2是640×480维的，图像分割时块大小取16×16，迭代阈值取0.8。

为了衡量融合后图像效果，采用标准差、熵、平均梯度和交叉熵来评价融合后的图像效果，其中，标准差值、熵值和平均梯度值越大，表示融合的图像效果越好，越清晰，而交叉熵越小，则表示融合的图像质量越好。表1展示了本发明与平均值方法对六幅多聚焦图像的实验结果，对标准差值、熵值、平均梯度、交叉熵做了对比，其实验数据如表1：

表1本发明与平均值方法图像融合后的各性能参数对比

从表1中可以看出，本发明的融合结果都优于平均值方法的融合结果，本发明的标准差、熵和平均梯度都比平均值方法的大，说明本发明的融合方法从源图像中获取的信息丰富，细节表现力强，更接近标准图像，本发明的交叉熵比平均值法的交叉熵值小，说明本发明方法与标准图像的差异越小，获取的信息越多，丢失的信息越少，这说明本发明能够很好点保持图像的边缘细节信息，融合效果明显优于平均值方法。

将本发明与现有平均值方法对两幅标准测试图像Clock1、Clock2进行图像融合，融合效果对比图如图2所示，其中，图2(a)是512×512维的Clock1的原始图像，图2(b)是512×512维的Clock2的原始图像，图2(c)是在本发明的融合结果图像，图2(d)是平均值方法的融合结果图像。从图2中可以看出，本发明图像融合结果比平均值方法的图像融合结果清晰。

将本发明与现有平均值方法对两幅标准测试图像Pepsi1、Pepsi2进行图像融合，融合效果对比图如图3所示，其中，图3(a)是512×512维的Pepsi1的原始图像，图3(b)是512×512维的Pepsi2的原始图像，图3(c)是在本发明的融合结果图像，图3(d)是平均值方法的融合结果图像。从图3中可以看出，在图像中的人物头部部分，本发明方法融合结果图像轮廓清晰，没有模糊感，而平均值方法的融合结果图像人物头部轮廓不够清晰，模糊感很明显，对于钟表部分，本发明方法融合结果图像明显比平均值方法融合结果图像在数字部分清晰。

将本发明与现有平均值方法对两幅标准测试图像Lab1、Lab2进行图像融合，融合效果对比图如图4所示，其中，图4(a)是512×512维的Lab1的原始图像，图4(b)是512×512维的Lab2的原始图像，图4(c)是在本发明的融合结果图像，图4(d)是平均值方法的融合结果图像。从图4中可以看出，本发明融合结果图像比平均值方法的融合结果图像清晰。

2)本发明与平均值方法对图像聚焦区域融合的仿真实验

本实验的主要目的是对比本发明方法和平均值方法对于待融合图像聚焦区域的融合效果，通过峰值性噪比PSNR值来衡量两种方法的融合效果，PSNR值越大，说明图像越清晰，融合效果越好。其实验数据采用两幅标准测试自然图像：Pepsi1和Pepsi2，均为512×512维的。

图5(a)是待融合图像Pepsi1，其中矩形框中的部分为聚焦Pepsi1的聚焦区域；图5(b)是待融合图像Pepsi2，其中矩形框中的部分为聚焦Pepsi2的聚焦区域；图5(c)为本发明方法对于Pepsi1的聚焦区域融合效果图；图5(d)为本发明方法对于Pepsi2的聚焦区域融合效果图；图5(e)为平均值方法对于Pepsi1的聚焦区域融合效果图；图5(f)为平均值方法对于Pepsi2的聚焦区域融合效果图；从视觉效果上看本发明方法视觉效果好，图像清晰，信息保持度高，而平均值方法图像有模糊感，边缘信息保持度不好。

为了进一步说明本发明的优势，采用PSNR来对图像Clock1、Clock2、Pepsi1、Pepsi2的聚焦区域融合后的图像效果进行评价，PSNR值是衡量图像质量效果的一个重要指标，PSNR值越高说明图像越接近原图像，失真越少。表2列出了图像Clock1与Clock2，Pepsi1和Pepsi2的聚焦区域的融合结果图与原图像的PSNR值，其实验数据如下所示：

表2聚焦区域融合结果图与原图的峰值性噪比值

从表2可以看出，本发明的PSNR明显高于平均值融合方法，说明本发明融合后图像的聚焦区域更接近原图像，融合效果更好。

综上所述，本发明基于压缩感知框架下图像融合的特点，设计出了一个傅里叶系数低频全采的变密度观测模型，结合局部傅立叶谐波系数融合的特点，提出了压缩感知框架下的一种多策略图像的融合方法，与传统的融合方法相比具有计算复杂度低，图像清晰度高，细节信息保持较好，融合效果好的特点。

Claims (3)

1.一种基于压缩感知框架的多策略图像融合方法，包括如下步骤：

(1)分别输入原图像A和原图像B，并将原图像A和原图像B分为大小为C×C的局部图像X1和X2，C×C取8×8或16×16；

(2)对局部图像X1进行傅里叶变换得到傅里叶系数矩阵y1，对局部图像X2进行傅里叶变换得到傅里叶系数矩阵y2；

(3)采用傅里叶系数低频全采的变密度观测模型，对傅里叶系数矩阵y1进行观测得到观测向量f1，对傅里叶系数矩阵y2进行观测得到观测向量f2；

(4)根据观测向量f1和观测向量f2计算两幅图像的谐波系数H₁、H₂以及频谱匹配度S：

$H_1=\sqrt{\frac{\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}{\left|{f1}_{i,j}\right|}^2}{\underset{m=1}{\overset{C}{\mathrm{\Σ}}}\underset{n=1}{\overset{C}{\mathrm{\Σ}}}\left|X1(m,n)\right|}}$

$H_2=\sqrt{\frac{\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}{\left|{f2}_{i,j}\right|}^2}{\underset{m=1}{\overset{C}{\mathrm{\Σ}}}\underset{n=1}{\overset{C}{\mathrm{\Σ}}}\left|X2(m,n)\right|}}$

$S=\frac{2\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}\left|{f1}_{i,j}\right|\×\left|{f2}_{i,j}\right|}{\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}{\left|{f1}_{i,j}\right|}^2+\underset{i=1}{\overset{C}{\mathrm{\Σ}}}\underset{j=1}{\overset{C}{\mathrm{\Σ}}}{\left|{f2}_{i,j}\right|}^2}$

其中：H₁为图像A的谐波系数，X1(m，n)为局部图像X1在横坐标为m、纵坐标为n处的像素值，m∈[1，C]，n∈[1，C]，C为图像分块大小，f1_i，j为观测向量f1在横坐标位置为i、纵坐标位置为j处的值，i∈[1，C]，j∈[1，C]，H₂为图像B的谐波系数，X2(m，n)为局部图像X2在横坐标为m、纵坐标为n处的像素值，f2_i，j为观测向量f2在横坐标位置为i、纵坐标位置为j处的值；

(5)选取一个阈值T∈[0.6，1]，计算最大加权系数W_max和最小加权系数W_min：

W_max＝0.5-0.5×(1-S)/(1-T)

W_min＝1-W_max

其中S为步骤(4)中得到的频谱匹配度；

(6)根据步骤(4)中所得到的谐波系数H₁、H₂和频谱匹配度S，以及步骤(5)中所得到的最大加权系数W_max和最小加权系数W_min，计算融合后的观测向量f：

若S≥T且H1≥H2，则融合后的观测向量f＝W_max×f1+W_min×f2；

若S≥T且H1＜H2，则融合后的观测向量f＝W_min×f1+W_max×f2；

若S＜T且H1≥H2，则融合后的观测向量f＝f1；

若S＜T且H1＜H2，则融合后的观测向量f＝f2；

其中T为阈值，T∈[0.6，1]，f1为步骤(3)中得到的图像A的观测向量，f2为步骤(3)中得到的图像B的观测向量；

(7)根据步骤(6)中所得到的融合后的观测向量f，进行基于Wavelet-Curvelet框架的Split Bregman重构算法迭代20次，最后得到所需要的融合图像。

2.根据权利要求1所述的压缩感知框架下的多策略图像融合方法，其中步骤(3)所述的按照傅里叶系数低频全采的变密度采样的模型，对傅里叶系数矩阵y1进行观测得到观测向量f1，按如下步骤进行：

(2a)设定采样模型是取值仅为0或1的矩阵，将值为1的点作为采样点，根据输入图像A的大小设定矩阵B：若输入图像A的大小为m×m，则设矩阵B的大小为m×m且取值全为0，将坐标值为(0.5×m，0.5×m)的点a作为圆心，以0.3×0.3×m为半径的圆，将圆内所有点都作为采样点，将矩阵B中这些位置上的点值设为1；

(2b)按照

的采样密度对步骤(2a)圆外的采样点随机选取，取到的点值设为1，未取到的点值设为0，其中r为圆外任取一点b到圆心a的距离长度，R是以圆心a为起点通过点b，以图像边缘为终点的线段长度；

(2c)重复步骤(2a)和步骤(2b)z次，取z次的平均值作为最终的观测模型M，其中z取值为500到1000；

(2d)用步骤(2c)中得到的观测模型M分别点乘步骤(2)中得到的傅里叶系数矩阵y1得到观测向量f1。

3.根据权利要求1所述的压缩感知框架下的多策略图像融合方法，其中步骤(3)所述的按照傅里叶系数低频全采的变密度采样的模型，对傅里叶系数矩阵y2进行观测得到观测向量f2，按如下步骤进行：

(3a)设定采样模型是取值仅为0或1的矩阵，将值为1的点作为采样点，根据输入图像A的大小设定矩阵B：若输入图像A的大小为m×m，则设矩阵B的大小为m×m且取值全为0，将坐标值为(0.5×m，0.5×m)的点a作为圆心，以0.3×0.3×m为半径的圆，将圆内所有点都作为采样点，将矩阵B中这些位置上的点值设为1；

(3b)按照

的采样密度对步骤(3a)圆外的采样点随机选取，取到的点值设为1，未取到的点值设为0，其中r为圆外任取一点b到圆心a的距离长度，R是以圆心a为起点通过点b，以图像边缘为终点的线段长度；

(3c)重复步骤(3a)和步骤(3b)z次，取z次的平均值作为最终的观测模型M，其中z取值为500到1000；

(3d)用步骤(3c)中得到的观测模型M分别点乘步骤(2)中得到的傅里叶系数矩阵y2，得到观测向量f2。

Images