一种基于多尺度空间正则化模型的超分辨率制图方法

摘要

本发明公开了一种基于多尺度空间正则化模型的超分辨率制图方法。该方法利用软分类所获得的组分图像作为输入数据，极小化组分残差，通过细尺度的空间匀质性作为光滑项，同时利用粗尺度的空间异质性进行惩罚，从而使得产生的超分辨率制图结果不仅有光滑的边缘而且能够较好的保持局部细节(如连续的地物类别)，不仅可以消除斑点问题，更重要的是可以很大程度的避免常见的过光滑现象。本发明具有实用性广，计算有效，特别对于小放大因子具有较优的超分辨率制图精度等优点，适用于土地覆被/土地利用制图、遥感影像目标提取及遥感影像中地物变化检测等利用遥感图像精细制图的相关领域。

说明书

技术领域

本发明属于空间信息技术，图像处理和信号处理的交叉领域，具体涉及一种基于多尺度 空间正则化模型的超分辨率制图方法。

背景技术

利用数字遥感影像制作专题地图(影像分类)是遥感技术研究中的热点领域。然而，由 于地表复杂性，光谱传输过程中各种环境的影响及传感器探测元件本身的物理特性等因素的 作用，遥感影像中的像元很少是由单一地物组成，一般都包含两种或两种以上的感兴趣地物。 这种由多种地物组成的像元被称为混合像元，其光谱特征是所对应的不同土地覆盖类型光谱 响应的综合特征。混合像元在栅格编码影像中是普遍存在的，所占比重也很大，尤其在不同 类型地物的交界地带最为普遍，它的存在使同一类型或属性地物的遥感影像灰度值具有一定 的变化范围。解决混合像元的问题能够进一步细化影像分类结果，使遥感技术应用进入到亚 像元尺度，更加详细地表示地物细节信息，而不仅仅是停留在像元尺度上，就成为一个必要 而且有意义的研究内容。

由于超分辨率制图逆问题的本质，基于正则化方法能有效获得一个近似最优的超分辨率 制图结果。在基于正则化框架的超分辨率制图方法中，按照待优化的目标函数可以分为两类： 一类是保持谱信息，最有代表的是Kasetkasem于2005年提出的基于马尔可夫模型的方法。 该方法以谱数据似然概率作为保真性，空间依赖性表示的先验概率作为正则项。还有Xu等 人于2013年提出的基于观察模型的拟合误差作为保真性，某种先验作为正则项。然而虽然 基于光谱数据保真性的方法有同时利用光谱和空间约束的优势，但由于集成了解混过程，解 混的不确定性必然影响超分辨率制图的结果，且使得制图结果的解释过于复杂。另一类是 Ling等人于2014年提出的保持组分信息，不受解混的影响。该方法极小化输入组分数据和 结果制图之间的差异，并利用细尺度的空间依赖性作为正则项来光滑化。

Ling等在2014年提出的保持组分信息并利用细尺度空间依赖性表示正则项的超分辨率 制图方法命名为SRM_REG，该方法优点在于：(1)仅利用组分信息，不需要辅助数据；(2) 该方法可以获得光滑的边界。但仍然有待改进：容易产生过光滑现象，使得局部细节被周围 大多数类别同质化而失真。针对该问题，本发明提出一种基于多尺度空间正则化的超分辨率 制图方法(SRM_MSSR)，延续了SRM_REG的正则化框架，并进一步提出利用亚像元在粗 尺度的空间异质性来惩罚由于仅考虑细尺度的空间依赖性所产生的过光滑问题。利用惩罚 项，SRM_MSSR能更好的保留细节，如细条状地物的连续性等。

发明内容

本发明要解决的技术问题为：提供一种惩罚由于细尺度空间相关性产生的过光滑问题的 超分辨率制图方法，该方法通过粗尺度的空间异质性对细尺度的光滑效应进行惩罚来抑制过 光滑问题。

本发明解决上述技术问题采用的技术方案为：一种基于多尺度空间正则化模型的超分辨 率制图方法，其特征在于：该方法包括如下步骤：

步骤1、预处理遥感数据，采集样本点，进行软分类获取被处理影像上每个像元的软分 类信息，即每个像元内部的各组分也称作端元的面积矢量或称丰度；

步骤2、设置超分辨率制图的放大系数N,N＝2,3,…,为自然数；

步骤3、基于空间吸引模型生成初始超分辨率制图结果；

步骤4、设置正则参数λ1和λ2，初始温度及最大迭代次数ite_max；

步骤5、基于多尺度空间正则化的优化模型对每个亚像元的类别进行优化，具体过程为： 当迭代次数小于ite_max时，对于每个亚像元，在其父像元中出现的地物类别中随机选一个 指派给该亚像元；对于该亚像元计算目标函数；估计该亚像元在随机指派前后的目标函数差， 如果目标函数减小，就接受该指派的类别为该亚像元的类别，否则利用模拟退火算法依据温 度参数以一定的概率，接受该指派；更新迭代次数和温度参数；

步骤6、重复步骤5至最大迭代次数；输出最终的超分辨率制图结果。

进一步的，所述步骤5中提出利用模拟退火算法以其获得近似全局最优解。

进一步的，所述步骤5中提出的优化模型既包括细尺度光滑项同时还包括在粗尺度抑制 过光滑的惩罚项，利用图像的组分信息计算粗尺度的空间异质性，从而惩罚由于仅利用细尺 度的空间依赖性所造成的过光滑问题，具体提出的优化模型形式化描述如下：

$X*=\arg \underset{X}{\min }E(X,D)+{\lambda }_1{R}_1\left(X\right)+{\lambda }_2{R}_2\left(X\right)$

$s.t.\underset{l=1}{\overset{c}{\Sigma }}{x}_{j,l}=1$

其中，X*表示待求解的最优超分辨率制图，表示超分辨率制图X同 输入的组分图像Y之间的差异性，表示细尺度亚像元的空间匀质 性，表示粗尺度的空间异质性，λ₁,λ₂为正则化参数，约束条件 表示每个亚像元只能取某一种地物类别。

本发明与现有技术相比的优点在于：

(1)、在仅利用组分数据的前提下，能够进一步提高遥感影像地物分类制图超分辨率重 构的精度，尤其是局部细节；

(2)、在仅利用组分数据的前提下，能够在很大程度上有效抑制遥感影像超分辨率地物 分类制图中仅考虑匀质性约束所造成的过光滑现象。

附图说明

图1为本发明一种基于多尺度空间正则化的超分辨率制图方法的实现流程图；

图2为本发明的实施过程所使用的模拟参考地物覆被制图；

图3为本发明的演示示例所使用的软分类结果；

图4为本发明的演示示例所使用的模拟初始超分辨率制图结果；

图5为本发明的演示示例所产生的中间制图结果；

图6为本发明的演示示例所产生的近似最优超分辨率制图结果。

具体实施方式

下面结合附图以及具体实施例进一步说明本发明。

如图1所示，一种基于多尺度空间正则化的超分辨率制图方法，在基于多尺度空间正则 化框架的超分辨率制图模型中，空间依赖性正则项包含两部分，一部分为光滑项，另一部分 为惩罚项。本发明的具体实施步骤如下：

步骤1，预处理遥感数据，采集样本点，进行软分类获取被处理影像上每个像元的软分 类信息，即每个像元内部的端元组分面积矢量或称丰度。假设在一个3×3像元的区域包含两 种土地覆被类型C1(白色栅格)和C2(灰色栅格)。通过软分类获得中心像元及邻域像元 的两种土地覆被类型所占的面积比例分别为Pm{0.25,0.75},P1{1,0},P2{1,0},P3{0.25,0.75}, P4{1,0},P5{0.75,0.25},P6{0.25,0.75},P7{0.75,0.25},P8{1,0}，如图3所示，其放大因子为2 的参考map如图2所示。

步骤2，设定放大因子为s＝2；

步骤3，基于空间吸引模型产生初始超分辨率制图结果，如图4所示。

步骤4，预设正则参数λ₁＝0.4和λ₂＝0.25，及初始温度T＝100和最大迭代次数ite_max＝30。

步骤5，顺序扫描中心像元Pm中的每一个亚像元p1-p4，并按公式(1)计算当前类别 的目标值及指派为其他类别时的目标值。

$L_j(X,d_i)=\underset{l=1}{\overset{c}{\Sigma }}{(y_{i,l}-d_{i.}·x_{.l})}^2+{\lambda }_1\underset{k_1\in N_j^{\left(1\right)}}{\Sigma }{\phi }_{j,k_1}^{\left(1\right)}\delta (j,k_1)+{\lambda }_2\underset{k_2\in N_i^{\left(2\right)}}{\Sigma }\underset{l=1}{\overset{c}{\Sigma }}{x}_{j,l}{\phi }_{j,k_2}^{\left(2\right)}{y}_{k_2,l})---\left(1\right)$

其中L_j(X,d_i)表示第j个亚像元的目标函数，X为超分辨率制图，d_i为第i个像元的 downsampling向量，c为地物类别个数(此例中为2)，y_i,l为第i个像元l类的丰度值，δ(j,k₁) 为第j个亚像元同邻域亚像元的kronecker delta函数，和分别为细尺度和粗尺度的邻 域权重，为第i个像元的粗邻域，为第j个亚像元的细邻域。

如果第j个亚像元的类别被修改为不同于当前的类别后L_j(X,d_i)减少，则更改该亚像元 类别为新类别；否则以与温度T相关的概率接受进行此次修改。

更新迭代次数，并按式(2)更新温度T：

T_ite＝T₀/log(ite+1)  (2)

重复步骤5，直到到达迭代次数。

举例说明公式(1)如何计算亚像元被指派不同类别后的目标函数值，此处用简化的邻 域权，即是邻域则权重为1，否则为0。

对于亚像元p1，当前类别为C2，在图4中，细尺度为C1类的邻域亚像元个数为6，粗 尺度邻域像元P1-P8的C2类组分和为0+0+0.75+0+0.25+0.75+0.25+0＝2，计算目标函数值：

L1＝0+0.4*6+0.25*2＝2.9

当亚像元p1的类别更改为C1，基于图4，其细尺度为C2类的邻域亚像元个数为2，粗 尺度邻域像元P1-P8的C1类组分和为1+1+0.25+1+0.75+0.25+0.75+1＝6，计算目标函数值：

L1＝[(0.25-0.5)²+(0.75-0.5)²]+0.4*2+0.25*6＝2.425

修改后的目标函数值变小，接受该修改，生成新的制图结果如图5。

类似，对于亚像元p2，当前类别为C2，在图5中计算目标函数值：

L1＝[(0.25-0.5)²+(0.75-0.5)²]+0.4*5+0.25*2＝2.625

当p2的类别更改为C1，基于图5计算目标函数值：

L1＝[(0.25-0.75)²+(0.75-0.25)²]+0.4*3+0.25*6＝3.2

修改后目标函数值变大，不接受该修改，保持当前制图结果如图5。

类似，对于亚像元p3，当前类别为C2，在图5中计算目标函数值：

L1＝[(0.25-0.5)²+(0.75-0.5)²]+0.4*5+0.25*2＝2.625

当p3的类别更改为C1，基于图5计算目标函数值：

L1＝[(0.25-0.75)²+(0.75-0.25)²]+0.4*3+0.25*6＝3.2

修改后目标函数值变大，不接受该修改，保持当前制图结果如图5。

类似，对于亚像元p4，当前类别为C1，在图5中计算目标函数值：

L1＝[(0.25-0.5)²+(0.75-0.5)²]+0.4*4+0.25*6＝3.225

当p4的类别更改为C2，基于图5计算目标函数值：

L1＝[(0.25-0.25)²+(0.75-0.75)²]+0.4*4+0.25*2＝2.1

修改后目标函数变小，接受该修改，生成新的制图结果如图6。

依次类推扫描完所有亚像元为止，更新迭代次数与温度T后再重复步骤5。

步骤6，输出近似最优超分辨率制图结果，如图6。

本发明未详细阐述部分属于本领域技术人员的公知技术。

以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员 来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等 同替换、改进等，均应包含在本发明的权利范围之内。