路志英 周庆霞 李 鑫

（天津大学电气自动化与信息工程学院，天津，300072）

引 言

太阳能是一种绿色环保的可再生能源，光伏系统可以实现对太阳能的有效利用，因此对光伏系统的预测研究逐渐引起了人们的关注[1]。然而，天空中云的生成、运动和消融会造成光伏系统输出功率不稳定[2]。因此，结合云图预测光伏系统超短期输出功率（或太阳辐射）的方法成为了当前主流的研究内容[3]。对于实时和超短期预测，结合地基云图进行预测是现在广泛使用的方法[4]。全天空地基云图是通过特定的拍摄设备对天空进行拍摄所获取的图像。由于全天空成像仪（Total sky imager,TSI）的设备结构特点导致采集的云图存在信息缺失，严重影响了云识别、云量计算及云量分布的描述。因此，对地基云图进行修复具有至关重要的意义[5]。

云图修复属于图像修复范畴，主要可分为两大类[6]：一类是针对小尺度缺损的基于结构的图像修复技术；另一类是针对大面积破损的基于样本块匹配的纹理合成算法。Bertalmio等[7]首次提出数字图像修复这一术语，该技术主要应用于破损油画或照片的修复、指定目标物的移除；在此基础上，Chan和Shen提出了全变分模型[8]和基于曲率驱动扩散模型[9]。此类修复方法在修复小尺度破损图像时有很好的修复效果，但对大面积纹理区域缺失的图像进行修复时会产生模糊现象，因此主要适用于非纹理图像和破损程度不严重的纹理类图像。Criminisi等[10]于2004年提出一种基于样本的图像修复算法，该算法在确保可以修复较大破损区域的同时，还一定程度上保留了图像的结构和纹理信息，使得对物体的移除和区域的填充具有良好的效果。后续大量的专家学者对Criminisi算法进行研究，提出了大量的改进算法[11-14]。例如，彭坤杨等[11]通过基于图像平均灰度值的方法，缩短算法执行时间，修复效果与Cri⁃minisi算法持平，但丢失了图像的彩色信息；朱霞等[12]采用基于颜色区域分割的方法，其分配块的搜索限定在源样本块所覆盖的颜色区域中，提高了算法性能，但增加了算法的复杂度；韩明珠等[13]通过引入调节因子减小置信度对优先权的影响，间接增大数据项对优先权的作用，在一定程度上增加了算法的复杂度。

地基云图纹理细节比较丰富，本文在Criminisi算法的基础上提出了一种基于改进Criminisi算法的地基云图修复方法，即根据地基云图独特的红蓝比特征，对优先权函数进行改进，使修复顺序得到了改善，同时对搜索匹配块时的搜索区域进行了调整，使修复结果更加自然，提高了修复的准确性。

1 Criminisi算法简介

Criminisi修复算法的原理是在最佳样本上通过寻找待修复块的最佳匹配块及复制来进行图像修复。该算法主要由计算待修复块优先级、寻找最佳匹配块并填充以及更新置信度3部分构成。令Criminisi修复算法中符号标记如图1所示。其中：I代表整个图像；Ω代表图像I中待修复的区域；δΩ代表待修复区域与源区域的边界；Φ代表源区域，可为待修复区域提供最佳匹配样本，且Φ=I-Ω；ψp为位于待修复区域边界δΩ上大小为(2ε+1)×(2ε+1)的待修复块，中心点为p，一般情况下选择M×M（M=9）的方形窗口模块。

1.1 计算待修复块的优先级

Criminisi修复的核心与重点是不同待修复块优先级的计算与比较，图像源区域的每个像素点都有其不同的置信度，而待修补块由于其内部像素点的差别而具有不同的优先级。

对于以边界线∂Ω上的点p为中心的待修复块ψp，其优先权P(p)的计算公式为

图1 Criminisi算法标记图Fig.1 Criminisialgorithm marker map

式中：C(p)表示置信度项，用来衡量待修复块ψp中可靠信息的数量；D(p)为数据项，表示当前点p所达到的等照度线的强弱。C(p)和D(p)可由式(2，3)计算获得。

需要对C(p)进行初始化处理，当∀p∈ Ω时，置C(p)=0；当∀p∈ Φ时，置C(p)=1，有

式中：np表示p点处垂直于边界线δΩ的法向量；∇Ip⊥表示点p的等照度线方向；α为归一化参数，一般选择α=255。

根据式(1)计算待修复区域边界上所有点的优先级，确定出优先级最大的点所在的待修复块。

1.2 寻找最佳匹配块并填充

寻找最佳匹配块是在源区域中搜索与具有最大优先级的待修复块纹理最相近且与之匹配的块；填充则是把搜索到的最佳匹配块中的像素对应复制到待修复块的未知像素点，从而使未知像素点转化为已知像素点。

不同的匹配块与具有最大优先级的待修复块之间的匹配度公式为

式中：Ψp表示待修复块；Ψpˆ表示匹配块，其中 pˆ表示匹配块中的像素点；d(Ψpˆ,Ψp)表示 Ψp与 Ψpˆ块区间的差距。基于像素灰度平方差之和（Sum of squared differences，SSD）最小的匹配准则，待修复块与匹配块内像素之间存在的颜色差距可描述为

式中函数R()、G()、B()分别表示M×M（本文设定M=9）图像块中点(i,j)的红、绿、蓝三原色。

确定最佳匹配块之后，将搜索到的最佳匹配块中的像素对应复制到待修复块的未知像素点，实现待修复块的修补过程。

1.3 更新置信度

待修复块的修复实现之后，其中的边界点转化为内部的已知点，原待修复区域内的点变为已知点或者边界点，此时需要更新已知像素点的置信度以及边界点的优先权。将这些点的置信度统一设为与p点相同的置信度

对以上3步进行循环，直至待修复区域完成全部修复。

2 基于TSI设备的图像采集

2.1 图像采集设备

图2 全天空成像仪TSI-880外观Fig.2 TSI-880 appearance

TSI-880是美国Yankee Environmental Systems公司研制的全自动、全彩色天空成像系统，能够满足各种天气条件下对局地天空的观测需要。图2所示为TSI-880的外观，该设备的工作原理是通过CCD镜头向下观测一个半球型旋镜，捕捉从天空反射到镜上的图像；在球面镜上方有一条遮光带随太阳的运动而转动，用于避免阳光直射入镜头，起到保护TSI的作用[15]。TSI-880采集的地基云图为640像素×480像素的JPEG图像，通过网络服务器实时显示。

图3为TSI-880输出的云图图像，中部圆形区域为有效区域，因此处理图像时首先将圆形区域外部置白，以消除外部区域对图像处理的影响；其次，中部圆形区域中的遮光带和镜头支臂影像为待修复区域，对其进行置黑，处理后的结果如图4所示。

图3 TSI-880输出的云图图像Fig.3 Cloud image outputted by TSI-880

图4 待修复区域的显示Fig.4 Displayed area to be inpainted

2.2 遮光带和镜头支臂影像定位

TSI设备的镜面上存在保护其摄像头不受太阳强光损坏的遮光带，且设备用于支撑镜头的机械臂也在镜面上形成投影，使获得的云图信息不完整，影响图像的分析，因此为了修复地基云图还原云图真实情况，需要确定出遮光带和镜头支臂影像在云图中的位置。

遮光带在云图中的位置会随着太阳位置的变化而变化，因此需要计算太阳光斑在地基云图中的位置以对遮光带进行定位，以此确定待修复区域。在天文学中，可通过太阳天顶角和方位角确定太阳的位置。遮光带与镜头支臂影像定位的计算过程如下：

（1）获取TSI云图拍摄的时间和设备所在位置的经纬度。

（2）计算太阳赤纬角δ为

式中N为一年中的日期顺序序列号，从一年中的第一天（1月1日）到最后一天（12月31号），取值范围是1～365。

（3）计算太阳时角ω为

式中t表示真太阳时，以24 h计。

（4）计算太阳方位角γs和天顶角θz为

式中φ表示全天空成像仪TSI所处位置的纬度。

（5）确定遮光带位置。由于遮光带区域的轴线和云图中心线的夹角与太阳方位角γs（如图5所示）一致，所以可以通过太阳方位角来确定遮光带区域的中心线在云图中的位置。

（6）镜头臂的位置在图像中固定，可直接确定。

3 Criminisi图像修复算法的改进

3.1 优先权的改进

Criminisi算法在确保可以修复较大破损区域的同时，还一定程度上保留了图像的结构和纹理信息，使得对物体的移除和区域的填充具有良好的效果。然而由于地基云图纹理细节比较丰富，具有独特的红蓝比特征，因此将该特征引入置信项中，可实现对待修复块优先级的计算。

图5 方位角在地基云图中的显示Fig.5 Display of azimuth in the ground-based cloud image

在优先权函数中的置信项中加入了地基云图独特的红蓝比特征，避免了置信度快速下降造成的影响，使得图像修复时的顺序更加合理。对Criminisi算法的优先权计算公式作如下改进

式中RB(p)表示改进的置信度项，则RB(p)为

式中：ψp表示待修复块，选取9×9的方形窗口模块；rb(q)表示q点的红蓝比，即像素点q的红色通道和蓝色通道的分量比值。

3.2 基于启发信息的最佳匹配块获取

在Criminisi算法中每个待修复块采用遍历的方法对整个源区域进行全局搜索，这种搜索方法会增加算法的时间复杂度，使效率降低。而一块图像与周围相邻图像的相关性与距离成反比，采用全局搜索有时会得到离待修复块较远且相关性较低的匹配块，即过大的搜索范围可能会导致最终得到的修复图像与其邻域图像有较大的差别。

通常一幅图像具有很强的空间冗余度，而且相邻像素值的变化小，这些具有相关性的点形成了图像的结构和纹理。对于待修复块可以只选取待修复点周围的源区域作为待匹配区域，这样可以缩短搜索时间，同时保证修复图像在其邻域图像内搜索。为此，本文提出基于启发信息选择匹配区域大小的方法，步骤如下：

（1）待修复区域为遮光带和镜头支臂影像区域；

为匹配区域，如图6所示。图中：p为待修复点，w为遮光带宽度，M为以点p为中心的待修复块的边长，L为以点p为中心的正方形的边长，图中阴影部分即为匹配区域；

图6 待修复点的匹配区域Fig.6 Matching area of the point to be inpainted

（3）在匹配区域中遍历，寻找最佳匹配块。

3.3 改进后算法流程

改进后的具体实现步骤如下（流程图如图7所示）：

（1）对提取出的遮光带区域和镜头臂区域置黑；

（2）提取置黑区域的边界，并判断边界是否为空；若是，跳到步骤（7）；否则，继续下一步；

（3）计算边界上各点的优先级，确定修补的先后顺序，找到优先级最大的像素点所在的待修复块，先对其进行修复；

（4）确定待匹配区域。以待修复点为中心选取边长为L的正方形，该区域与源区域相交的部分即为匹配区域；

（5）计算待修复块和最佳匹配块的SSD，选取SSD最小的匹配块对待修复块进行修复；

（6）更新破损区域，返回步骤（2）；

（7）完成修复，输出修复图像。

4 实验结果与分析

图7 算法流程图Fig.7 Algorithm flowchart

应用本文算法对地基云图进行了修复，并与Cri⁃minisi算法、文献[13]算法的结果进行了比较。上述算法均在WIN 7 platform with Intel®Core TM i3-2120 CPU@3.30 GHz的配置条件下，利用Matlab 2010编程完成的。

图8 采集图像1Fig.8 Capture image 1

图8 是拍摄时间为2016-02-17-151930的地基云图。从图8(c,d)不难看出：虚线框中的遮光带区域的修复结果出现了块匹配错误，产生断层现象，造成结构不连续、不自然；而图8(e)相应区域修复后的图像较连续，没有明显的块匹配错误，保证了合理的填充顺序，达到了较好的视觉效果，而且自然可信，表明了其修复效果明显优于原始算法。

图9、图10为两组不同时刻地基云图的修复实验，在Criminisi修复算法及文献[13]修复算法得到的实验结果中可以看出结构不连续不自然的问题，而改进后的Criminisi修复算法得到的实验结果可以取得更好的效果。

图9 采集图像2Fig.9 Capture image 2

图10 采集图像3Fig.10 Capture image 3

表1给出了修复不同时刻的地基云图所用的时间。从表中数据可见：改进后的Criminisi修复算法的修复速度相对于Cri⁃minisi算法提升了5倍左右，文献[13]修复算法引入调节因子使得需要的修复时间更长，说明改进后的Criminisi修复算法降低了时间复杂度，提高了修复效率。

为了更好地反映改进算法的有效性，本文采用峰值信噪比（Peak signal to noise ratio,PSNR）评价指标对修复后的地基云图进行了客观评价，同时也对本文的比较方法进行了评价。PSNR计算公式为

表1 算法执行时间比较Tab.1 Algorithm execution time comparison

本文随机选取了图10（a）的局部区域进行置黑，并用本文算法对其进行修复，实验效果如图11所示，实验结果如表2所示。表中PSNR用于衡量修复结果相对于原图的保真度。一般来说PSNR值越大，修复效果和质量越好，图像失真越小。综上所述，从PSNR值可以看出改进后的算法PSNR值提高；从修复的视觉效果上可以看出：改进后的算法比原始算法更加自然。

表2 不同算法修复结果的PSNR值Tab.2 PSNR values of differ ent r epair ed results

图11 局部区域修复效果Fig.11 Inpainted effect of local area

5 结束语

本文针对TSI拍摄的地基云图中存在的遮挡问题，结合地基云图特征，对Criminisi算法进行了改进，提出了一种基于改进Criminisi算法的地基云图修复方法，具体包括：

（1）引入地基云图独特的红蓝比特征作为置信项，将图像更多的已知信息融入到优先权计算中，使得含有更多信息的像素块具有更高的优先级，解决了传统Criminisi算法在计算优先权时由于置信度急剧下降为零所造成的修复顺序发生错误的问题。

（2）提出基于启发信息选择匹配区域大小的方法，避免了搜索到离待修复块较远的相关性较低的匹配块，降低了算法的时间复杂度。

实验表明，与传统Criminisi算法相比，本文修复算法可以在较短的时间内实现对不同时刻地基云图的修复工作，并取得更为满意的修复效果，为建立光伏发电预测模型奠定了基础。