2012年1月9日，我国成功发射第一颗自主的民用高分辨率光学传输型立体测绘卫星—资源三号(ZY-3) 卫星。2016年5月30日11时17分，我国成功发射资源三号02星，它是资源三号系列卫星的第2颗，实现两颗资源三号测绘卫星组网运行。资源三号系列卫星影像在空间分辨率、定位精度与时效性等方面代表了我国自主民用遥感卫星的领先水平，可以服务于基础测绘、国土、农业、环境、减灾、规划等各行业影像数据需求，具有广阔的应用前景。

然而对于绝大多数基于光学传感器的遥感影像来说，最大的问题是易受恶劣天气的影响，当云、雾或者水汽等在信号传播路线上形成一定规模，会导致图像的不准确甚至被遮盖，导致地物信息缺失，阻碍遥感影像发挥应有的作用价值，造成巨大的经济损失。例如2006、2007年北京全城有云覆盖的天气约占总量的46%，我国地处亚热带的南方地区云覆盖频率更高。通常，云雾覆盖率超过10~15%的影像不能应用于测绘研究中，然而根据1991年至1998年的卫星资料显示，由于云雾天气的影响，我国南方地区获取的NOAA/AVHRR气象卫星数据的平均有效利用率不足7%^[1]。由此可见，去云不仅是遥感影像进行准确解译的基础，也是增强遥感数据有效性的重要途径，故遥感影像去云具有十分重要的实际意义。

目前，遥感影像去云主要分为薄云去除和厚云去除，薄云去除的方法主要有基于直方图匹配去云法^[2-3]、基于多光谱图像去云法^[4-5]和基于单幅图像去云法等，厚云去除的方法主要有基于支持向量机法^[6-7]、基于光谱及几何法^[8]、基于数据融合去云法^[9-10]等。在遥感影像去云方法中，去除厚云的处理较去除薄云的方法更为棘手，引起了不少学者的关注。近几年，在压缩感知理论发展的基础上，利用过完备字典进行稀疏表示的方法在图像处理邻域得到了广泛的应用，也为遥感影像云去除处理提供了新思路。文献[11]提出的K-SVD字典训练算法，依据误差最小原则，对误差项进行SVD分解，选择使误差最小的分解项作为更新的字典原子和对应的原子系数，经过不断的迭代从而得到优化的解，实现输入信号字典自适应的稀疏表示。

本文应用多时域资源三号遥感影像，以K-SVD字典学习算法为基础，为其字典训练过程添加排序准则，使得不同图像的字典原子排列顺序尽可能相似，减少同区不同时的遥感影像间的差异，达到较好的去除云与阴影的效果。此方法对于薄云、厚云均适用，但针对厚云效果较为明显，故以下将以厚云情况为例。

1 K-SVD理论基础 1.1 压缩感知

压缩感知(compressed sensing, CS)也称为稀疏采样，是一个全新的采样理论，它通过开发信号的稀疏特性，在远小于Nyquist采样率的条件下，用随机采样获取信号的离散样本，然后通过非线性重建算法完美地重建信号^[12]。压缩感知理论已在信息论、图像处理、地球科学、光学、微波成像、模式识别、地质等领域得到广泛应用。

压缩感知理论的实质是利用多个观测值无失真地恢复出长度为M的原始信号，主要内容涉及信号的稀疏表示、测量矩阵的设计和信号恢复算法，其数学表达为：在已知长度M的一维测量值y和测量矩阵H(M×N)的基础上，求解欠定方程组y=Hx，得到长度N、稀疏度K(即含有K个非零值)的一维原信号x。但是一般的自然信号x本身并不稀疏，需要某种转换进行稀疏表示，故令x=Sα，S(N×L，L≫N)为稀疏基矩阵，α为稀疏系数(α只有K个是非零值(K≪N))，原理表示如图 1所示。

1.2 K-SVD字典训练

字典训练以压缩感知理论为基础，是数据稀疏表示的最优表示，能够训练出性能优良、规模更加紧凑的新字典，使字典的原子特性更加符合需要表示的图像信号。原来的关注点一直聚集在给定一个字典条件下分解信号的追踪算法，这个字典要根据一个适应一系列训练信号的特定的线性变换得来。而文献[11]证明了用于稀疏表示的过完备字典不仅可以根据给定的函数集选择得到，还可以依据输入的信号自适应的表示，此外，K-SVD具有很强的适应性，可以与任意一种常见的稀疏编码追踪算法结合，例如，基追踪(basic pursuit，BP)算法^[13]、匹配追踪(matching pursuit，MP)算法^[14]、正交匹配追踪(orthogonal matching pursuit，OMP)算法^[15-16]和焦点欠定系统求解^[17-18](the focal underdetermined system solver，FOCUSS)等。K-SVD的本质是使用包含原信号原子的超完备字典，由原子的稀疏线性组合恢复信号。它是由K-means聚类过程产生，可以基于现有字典进行稀疏编码，然后通过循环迭代逐个更新字典的原子和稀疏系数，字典列的更新与稀疏表示的更新能使最终输出的字典更好地适应数据，并且聚合更快。由此可见，K-SVD算法的目的是寻找最佳过完备字典表示W个样本信号{y_i}_i=1^W，它是稀疏表示的重要成分，求解满足式

式中，D=[d₁，d₂，…，d_L]∈R^M×L、A=[α₁，α₂，…，α_W]∈R^L×W分别代表所求的过完备字典和信号的稀疏矩阵；Y=[y₁，y₂，y₃，…，y_W]∈R^M×W为已知的数据集；T为给定的稀疏系数中非零元素个数的上限值。

对于尺寸为V×V的图像，可以被滑动重叠的分为若干尺寸为b²，b≪V的小块。每一个小块都可以表示为b²维向量y∈R^M×1，该向量就是原始图像的一个块向量。字典D是一个M×L的矩阵，D∈R^M×L，由L个原子组成，每个原子同样也为M维向量d_i∈R^M×1，由于字典D是过完备的基的集合，所以M≪L。字典学习和稀疏表示都是对于这些小块进行的，决定了字典的尺寸，使得字典的行数为M=b²，以及字典的列数(训练的原子个数)由冗余因子r和b共同决定为L=rb²。按照K-SVD通常设定，本文试验参数取r=4，b=8，得到块大小b²=64，字典行数M=64，原子个数L=256，K-SVD迭代次数取10，最大训练块的个数设置大于(V-b+1)²，并且根据测试取稀疏表示全局误差参数为25×1.15，基于块恢复的流程图如图 2所示。

该算法运行的主要步骤是：① 初始字典：设置以标准正交基为列的初始字典D⁰(多为DCT矩阵)；② 稀疏编码：选用任意一个追踪算法计算每一个样本信号y_i对应的稀疏向量α_i，其求解满足(i=1，2，…，W)，得到新的稀疏向量α_i^*；③ 更新字典原子：利用原子d_k的全局表示误差矩阵E_k=Y-，并应用奇异值分解(singular value decomposition，SVD)原子的误差矩阵，即E_k^*=UΔV^T，则U的第1列为更新的字典原子d_k^*, V的第1列乘以Δ(1，1) 即为更新的系数向量α_k^*；④ 循环步骤② 和③，直至满足误差阈值。

近年来，K-SVD算法已在图像的压缩、识别、修饰、去噪和重建等领域取得了突出效果，并被不断优化改进，但在遥感影像去云方面刚开始起步，为遥感影像去云提供了新的思路。文献[19]应用压缩感知理论重建云区域，分别试验比较了不同尺寸的云覆盖和不同类型地区的云区域恢复效果；文献[20]首次提出将定量遥感影像与字典学习结合；文献[21]将多光谱遥感影像结合多时域字典学习重建云污染区域。本文在改进K-SVD字典学习的原子训练过程的基础上结合正交匹配追踪算法，重建遥感影像云区域，丰富地物信息。

2 本文算法 2.1 AO-DL算法

将字典学习算法应用于遥感影像云去除，用字典稀疏表示同区不同时的遥感影像时，由于不同的拍摄条件和部分不同的图像内容(地物的变化、云的遮盖)，得到的字典原子的排序绝对不同，因此用有云影像的字典和无云影像的稀疏系数矩阵直接恢复有云图像难以实现。但是若对两个字典的原子按照特定的规则重新排序，使其顺序尽可能相同时，则可以减少影像间的差异，能很大程度地改善影像去云的效果。

本文改进了字典学习算法，在基于以上块训练所得字典的基础上，再将字典原子按照某种特定的准则进行重排序，得到基于块的字典训练原子重排序，称为变序字典学习(atoms-reordered dictionary learning，AO-DL)。目前，字典学习算法有Sparsenet字典学习算法^[22]、MOD(method of optimal directions)字典学习算法^[23]、K-SVD字典学习算法、在线字典学习(online dictionary learning，ODL)算法^[24]和结构字典学习算法^[25]等。本文选用K-SVD算法和以最大相关系数为重排序准则实现AO-DL，其中K-SVD算法已在1.2节中简单介绍，以下将重点说明字典原子重排序的方法。

应用最大相关系数进行字典原子重排序：

(1) 数据准备。寻找一张灰度均匀的图像作为标准图，记为image0；输入图，记为image(i)。

(2) K-SVD计算，得到图像对应的字典。image0对应的K-SVD图像字典记为D0_M×L，image(i)对应的K-SVD图像字典记为D(i)_M×L。

(3) 关联矩阵计算。依次计算字典D0的L个原子对应字典D(i)的每个原子的相关系数(correlation coefficient)，形成系数矩阵，记为CM_L×L。即字典D0的第n列原子与字典D(i)的所有原子的L个相关系数对应为CM第n列的L个元素。

(4) 计算排列顺序，记为order。根据CM每列中最大相关系数所在的位置，形成一个排列顺序。也就是说，order是一个1×L的数组，它的第n个元素即为CM第n列的最大相关系数所在的那个行数f，即order[n]=f。这就意味着字典D(i)的第f列原子与D0的第n列原子是最相关的。

(5) 得到原子重新排序后的字典D(i)_M×L^*。按照order元素顺序对字典D(i)的原子进行重新排序，得到新字典D(i)_M×L^*。

由此，输入图像的字典基于标准图像的字典按照一种的有规则的排列顺序(本文所指的规则就是最大相关法)生成了新字典，具体AO-DL算法实现流程图如图 3所示。

2.2 AO-DL云去除

AO-DL算法完成后，可以令输出的字典按照一定的规则进行原子重排序。将此算法应用于遥感影像云去除中，输入一张无云影像作为参考图像img_A和一张有云影像作为目标图像img_B，两张图像经过AO-DL算法后对应产生两个带有各自影像属性以及相同的原子排序规则的字典，分别为字典D_A^*和字典D_B^*。然后，运用正交匹配追踪(OMP)算法令img_A中的部分无云区域图用字典D_A^*原子的线性组合稀疏表示出来，得到稀疏矩阵α。最后，将带有目标图像属性的字典D_B^*与稀疏矩阵α相乘，得到目标图img_B中对应有云区域的去云重建图img_b^*，算法流程如图 4所示。应用该方法得到的去云图可以在拥有原影像属性的同时有效地去除云和阴影的遮盖。

3 试验与结果 3.1 数据准备

论文试验数据列举两组资源三号(ZY-3) 卫星影像，第1组为分别于2015年2月9日与2015年8月15日拍摄的天津市宁河区以农田区域为主的影像，取4763×3576像素，1个像素长2.1 m；第2组为分别于2015年8月15日与2016年2月8日拍摄的天津市滨海新区影像，取3046×2500像素，1个像素长2.1 m。试验数据集如图 5所示，详细信息见表 1。

3.2 试验

论文主要讨论遥感影像云去除处理，对于云区域的选取采用手工框选实现的，因此本文不涉及云覆盖区域的自动监测内容。在算法试验之前，应用ENVI先对同区域不同时域的两幅原始图像进行影像配准，满足RMS不超过1.0，然后对两幅图像的重叠区域进行矩形剪裁。由于本文算法只需要保证框选出的试验矩形图中包含云即可，而不用像一般检测算法一样精确地检测出云的具体轮廓。本文试验矩形区域的选取采用以相同的经纬度坐标矩形角点为参照标准，固定长宽，框选出试验矩形范围，即有云图像与无云图像的试验矩形的角点经纬度坐标是同一的，从而保证试验矩形区域像素的对应性。试验从原始有云污染的遥感影像中截取云覆盖区域以及无云影像的对应区域，并按照图 4所示流程进行试验, V组数据K-SVD算法的重建效果如图 6所示，两组数据云覆盖区域本文改进算法AO-DL的重建效果如图 7和图 8所示。由试验结果可见，变序字典学习AO-DL算法对于资源三号遥感影像去云有着较好的效果。

采用均方根误差(root mean square error，RMSE)定量描述重建图像的恢复效果，计算有云影像的局部无云区域与重建图像的对应区域的误差，例如图 7效果与图 6效果比较，AO-DL方法重建图像局部区域与对应原清晰影像的均方根误差为2.286，而未改进的字典学习法直接参与去云的均方根误差为4.524，可见将字典原子按照特定的规则重新排序后，重建图像的质量得到了显著的提升。改进后的方法与现有的绝大多数去云方法相比，其优点为：降低了云检测要求、操作简单、字典中含有原图像的属性等，其缺点为：算法处理时间较长，并且精度有待提高。

4 总结与展望

变序字典学习AO-DL算法应用多时域遥感影像，可以利用含有原图像属性的字典进行遥感影像的云污染区域重建，得到较好的恢复效果，并适用于不同类型区域。论文方法与现有绝大多数去云方法相比，有着降低了云检测要求、操作简单、字典中含有原图像的属性等优点，但是该方法在大面积成片云区域的情况下处理时间较长，且精度有待提高，主要因为该算法是在K-SVD算法的基础上进行的改进，应用K-SVD字典学习算法重建影像存在较原始清晰影像稍有模糊的问题，但对于有严重信息丢失的影像而言，此算法可以成功的补充大部分丢失的信息，因此相信字典学习算法将会更加广而应用。

后续主要有以下两点改进方向：一是进一步改善算法，减缓边缘模糊，提高精度；二是将尝试结合不同源的遥感影像加入到字典学习算法中，补充影像信息，从而提高云去除重建图像的精度。