胡 中，倪杭飛，李 偉，宓霄凌，白 帆，王伊娜

(浙江中控太陽(yáng)能技術(shù)有限公司，杭州 310053)

0 引言

在塔式太陽(yáng)能熱發(fā)電站中，由定日鏡鏡場(chǎng)將太陽(yáng)輻射反射到吸熱器上，然后吸熱器將接收的太陽(yáng)輻射能轉(zhuǎn)化成熱能，為后續(xù)的太陽(yáng)能熱發(fā)電提供能量。在該過(guò)程中，云層對(duì)定日鏡場(chǎng)上方太陽(yáng)輻射的遮擋是影響吸熱器接收的太陽(yáng)直接輻射(DNI)的主要因素，當(dāng)云層遮擋定日鏡上方的太陽(yáng)輻射時(shí)，定日鏡鏡場(chǎng)反射的太陽(yáng)輻射減少，吸熱器接收的DNI也相應(yīng)下降；而當(dāng)云層離開(kāi)時(shí)，定日鏡鏡場(chǎng)反射的太陽(yáng)輻射增大，吸熱器接收的DNI也隨之升高。因此，云層短時(shí)間的遮擋和離開(kāi)會(huì)使吸熱器接收的DNI發(fā)生隨機(jī)波動(dòng)，造成吸熱器表面的受熱情況發(fā)生劇烈變化。而通過(guò)超短期(0～30 min)DNI預(yù)測(cè)可以提前對(duì)定日鏡鏡場(chǎng)進(jìn)行調(diào)控，從而解決這一問(wèn)題。

基于上述原因，本文提出了一種基于極值梯度提升(eXtreme Gradient Boosting，XGB)算法及全天空成像儀圖像的超短期DNI預(yù)測(cè)方法。通過(guò)使用基于太陽(yáng)位置的全天空成像儀標(biāo)定方法來(lái)完成全天空成像儀的姿態(tài)校正和畸變矯正；然后對(duì)晴空背景參數(shù)進(jìn)行擬合并建立晴空庫(kù)，為云層的檢測(cè)提供基礎(chǔ)；根據(jù)云的運(yùn)動(dòng)矢量，利用外推法得到未來(lái)云層可能會(huì)遮擋太陽(yáng)的圖像點(diǎn)，并提取圖像特征，用于基于XGB算法的DNI超短期預(yù)測(cè)模型的訓(xùn)練，最后利用訓(xùn)練后的模型完成超短期DNI預(yù)測(cè)。

1 全天空成像儀

本文選用型號(hào)為EKO SRF-02的全天空成像儀，其分辨率為2274×1702，視野大于180°，可實(shí)現(xiàn)全天空監(jiān)控，其實(shí)物圖及安裝效果圖分別如圖1、圖2所示。將該全天空成像儀設(shè)置為每40 s采集2幀圖像，其中的一幀圖像為自動(dòng)曝光，另一幀圖像為欠曝光。本文后續(xù)均基于欠曝光圖像進(jìn)行處理。

圖1 型號(hào)為EKO SRF-02的全天空成像儀實(shí)物圖Fig.1 Physical image of all-sky imager model EKO SRF-02

圖2 全天空成像儀的安裝效果圖Fig.2 Installation effect image of all-sky imager

全天空成像儀在安裝時(shí)需進(jìn)行水平調(diào)準(zhǔn)，其水平傾斜度不能超過(guò)0.5°；由于全天空成像儀使用的是魚眼鏡頭，會(huì)存在畸變，因此需要進(jìn)行標(biāo)定。

本文開(kāi)發(fā)了一種基于太陽(yáng)位置的全天空成像儀標(biāo)定方法。該方法不需要通過(guò)標(biāo)定板對(duì)全天空成像儀進(jìn)行標(biāo)定，而是通過(guò)采集一整天的晴空?qǐng)D像，根據(jù)圖像上太陽(yáng)的中心位置，并結(jié)合理論太陽(yáng)高度角及太陽(yáng)方位角計(jì)算全天空成像儀的姿態(tài)參數(shù)及畸變參數(shù)，從而實(shí)現(xiàn)全天空成像儀的姿態(tài)校正與畸變矯正，減少了現(xiàn)場(chǎng)工作量。

2 建立晴空庫(kù)與云層的檢測(cè)

2.1 建立晴空庫(kù)

在對(duì)云層進(jìn)行檢測(cè)之前，需要建立晴空庫(kù)。最好在夏季采集晴空?qǐng)D像，因?yàn)橄募久刻斓奶?yáng)高度角的最大值接近一年中太陽(yáng)高度角的最大值；同時(shí)，需選擇大氣透射率最高的天氣。

同基于太陽(yáng)位置的全天空成像儀標(biāo)定方法一樣，建立晴空庫(kù)時(shí)也需要一整天的晴空?qǐng)D像，然后分別對(duì)每一幀圖像中的紅(R)、綠(G)、藍(lán)(B)通道進(jìn)行背景擬合，從而得到背景擬合參數(shù)。該背景擬合參數(shù)對(duì)應(yīng)該幀圖像的太陽(yáng)高度角，即每幀圖像的太陽(yáng)高度角對(duì)應(yīng)一組擬合參數(shù)，據(jù)此建立晴空庫(kù)。在進(jìn)行云層分割時(shí)，可以根據(jù)當(dāng)前時(shí)刻的實(shí)際太陽(yáng)高度角找到對(duì)應(yīng)的擬合參數(shù)(選擇與晴空庫(kù)中太陽(yáng)高度角差值最小時(shí)對(duì)應(yīng)的擬合參數(shù))，并根據(jù)該擬合參數(shù)生成晴空背景圖像。

根據(jù)文獻(xiàn)[1]，將晴空下天空亮度L建模為關(guān)于PZA(像點(diǎn)/天頂角，即像點(diǎn)對(duì)應(yīng)視線和天頂方向的夾角)和SPA(太陽(yáng)/像點(diǎn)角，即像點(diǎn)對(duì)應(yīng)視線和太陽(yáng)方向的夾角)的函數(shù)，具體可表示為：

同時(shí)，根據(jù)文獻(xiàn)[1]中基于實(shí)測(cè)數(shù)據(jù)的統(tǒng)計(jì)分析及建立的天空亮度仿真模型，可將晴空背景模型(下文簡(jiǎn)稱為“晴空模型”)擬合公式表示為：

式中：Ip(PZA,SPA)為天空亮度分布函數(shù)；φ1(PZA)為天空亮度漸變函數(shù)，簡(jiǎn)化記作F1；φ2(SPA)為散射指標(biāo)函數(shù)，簡(jiǎn)化記作F2。

其中：

式中：a1、a2、a3均為天空亮度漸變函數(shù)的擬合參數(shù)。

式中：a4、a5、a6、a7均為散射指標(biāo)函數(shù)的擬合參數(shù)。

在利用式(2)對(duì)全天空成像儀圖像進(jìn)行晴空背景擬合的過(guò)程中，發(fā)現(xiàn)該公式在部分情況下無(wú)法對(duì)太陽(yáng)周邊像素進(jìn)行精確擬合。

將式(2)兩側(cè)各減去φ1(PZA)，可得到：

由式(5)可知，Ip(PZA,SPA)–F1與F1呈線性關(guān)系。

B通道的Fg_b、Bg_b、F1及Fg_b–F1生成的曲線如圖3所示。其中：Fg_b為前景圖像B通道的值；Bg_b為背景圖像B通道的值。

圖3 B通道的Fg_b、Bg_b、F1及Fg_b–F1生成的曲線Fig.3 Curves generated by Fg_b，Bg_b，F(xiàn)1 and Fg_b–F1 of channel B

從圖3中可以看出，F(xiàn)g_b–F1與F1并非呈線性關(guān)系。這與根據(jù)式(5)得出的結(jié)論不相符，由此可知，式(5)在本文中并不適用。因此，本文參考式(1)、式(2)并結(jié)合圖3，修改了晴空模型，假設(shè)其是φ1(PZA)和φ2(SPA)的和，并認(rèn)為這2個(gè)函數(shù)的值獨(dú)立且無(wú)關(guān)聯(lián)，從而得到了新的擬合公式，即：

新擬合公式中用加法替代了式(2)(舊擬合公式)中的乘法，大幅降低了擬合計(jì)算過(guò)程的復(fù)雜度，使晴空模型得到了簡(jiǎn)化。

其中：

式中：b1～b3均為新擬合公式的天空亮度漸變函數(shù)的擬合參數(shù)。

式中：b4～b8均為新擬合公式的散射指標(biāo)函數(shù)的擬合參數(shù)。

利用新的晴空模型對(duì)晴空時(shí)的全天空成像儀圖像進(jìn)行擬合，得到擬合后的晴空背景圖像，如圖4所示，以及對(duì)應(yīng)的R通道的和的擬合結(jié)果，如圖5所示。

圖4 全天空成像儀圖像與擬合得到的晴空背景圖像Fig.4 Image of all-sky imager and processed clear sky background

圖5 利用新擬合公式得到的R通道的F′1、F′2值的擬合結(jié)果Fig.5 Fitting result of F′1 and F′2 values of R channel obtained by new fitting formula

對(duì)比了分別采用新、舊擬合公式得到的R、G、B通道的F′2的擬合結(jié)果，具體如圖6～圖8所示。

圖6 通過(guò)新、舊擬合公式得到的R通道的F′2值和F2值的擬合結(jié)果Fig.6 Fitting results of F′2 values and F2 values of R channel obtained by new and old fitting formulas

圖7 通過(guò)新、舊擬合公式得到的G通道的F′2值和F2值的擬合結(jié)果Fig.7 Fitting results of F′2 values and F2 values of G channel obtained by new and old fitting formulas

圖8 通過(guò)新、舊擬合公式得到的B通道的F′2值和F2值的擬合結(jié)果Fig.8 Fitting results of F′2 values and F2 values of B channel obtained by new and old fitting formulas

由圖6～圖8可知，針對(duì)本文所述型號(hào)的全天空成像儀拍攝的圖像，在太陽(yáng)高度角較小時(shí)，利用新擬合公式可以得到更好的擬合結(jié)果。

2.2 云層的檢測(cè)

云層的檢測(cè)范圍為天頂角80°以內(nèi)的區(qū)域。檢測(cè)前需要對(duì)該區(qū)域內(nèi)的地面建筑進(jìn)行裁剪并進(jìn)行插值修復(fù)，這是因?yàn)樵搮^(qū)域接近地平面，鏡頭畸變最為嚴(yán)重，沿圓弧方向進(jìn)行插值修復(fù)可以得到更好的插值結(jié)果。然后利用圖像變換公式分別計(jì)算前景NRBR(即歸一化的RB比值，此處為便于分析將值域映射到[0, 1])圖像和晴空背景NRBR圖像，并計(jì)算二者的差值圖像，即ΔNRBR圖像。NRBR的計(jì)算式可表示為：

由于晴空背景亮度的不均勻性，本文對(duì)ΔNRBR圖像進(jìn)行了歸一化。

在ΔNRBR圖像上的太陽(yáng)及其臨近區(qū)域會(huì)形成空洞，如圖9a所示，因此本文對(duì)該區(qū)域進(jìn)行了插值修復(fù)，插值范圍選擇SPA＜30°的區(qū)域，沿圓弧方向進(jìn)行線性插值，修復(fù)后的圖像如圖9b所示，其可用于準(zhǔn)確分割太陽(yáng)及其臨近區(qū)域的云層。

圖9 未處理與處理后的ΔNRBR圖像Fig.9 Images of unprocessed and processed ΔNRBR

使用固定閾值對(duì)歸一化的差值圖像進(jìn)行二值化，將較厚的云層(下文簡(jiǎn)稱為“厚云”)的分割閾值設(shè)為0.85，薄云的分割閾值設(shè)為0.5，此處的分割閾值為差值圖像灰度映射到[0, 1]后的取值。云的分割結(jié)果如圖10所示。

圖10 云分割前、后的圖像Fig.10 Images of before and after cloud segmentation

3 云的運(yùn)動(dòng)矢量計(jì)算

在對(duì)DNI進(jìn)行預(yù)測(cè)前需要先對(duì)云的運(yùn)動(dòng)矢量進(jìn)行計(jì)算。云的運(yùn)動(dòng)矢量是指云的整體運(yùn)動(dòng)方向和速度。本文使用歸一化相關(guān)系數(shù)匹配算法對(duì)全天空成像儀拍攝的前、后2幀圖像進(jìn)行匹配(匹配前需先對(duì)圖像進(jìn)行矯正，且只對(duì)天頂角80°范圍內(nèi)的像素進(jìn)行匹配，因?yàn)樵搮^(qū)域的像素失真較小)，并計(jì)算云的運(yùn)動(dòng)矢量，計(jì)算結(jié)果如圖11所示。圖中：綠線為相關(guān)系數(shù)匹配結(jié)果；紅線為計(jì)算得到的云的運(yùn)動(dòng)矢量，其中，線段長(zhǎng)度代表運(yùn)動(dòng)速度，箭頭代表運(yùn)動(dòng)方向。

本文使用卡爾曼濾波對(duì)不同時(shí)刻計(jì)算得到的云的運(yùn)動(dòng)矢量進(jìn)行了濾波，使計(jì)算結(jié)果更為穩(wěn)定。

圖11 云的運(yùn)動(dòng)矢量計(jì)算結(jié)果Fig.11 Calculation result of cloud motion vector

由于云的運(yùn)動(dòng)矢量是后續(xù)進(jìn)行云層分類及超短期DNI預(yù)測(cè)特征提取的基礎(chǔ)，因此，根據(jù)云的運(yùn)動(dòng)矢量，利用外推法計(jì)算未來(lái)30 min內(nèi)(時(shí)間間隔取1 min)云層可能會(huì)遮擋太陽(yáng)的30個(gè)圖像點(diǎn)。從上述30個(gè)圖像點(diǎn)提取圖像特征，進(jìn)行超短期DNI預(yù)測(cè)模型的訓(xùn)練與DNI實(shí)時(shí)預(yù)測(cè)。

本文的假設(shè)是：在短時(shí)間(0～30 min)內(nèi)云的形態(tài)保持不變，且云的運(yùn)動(dòng)形式為整體平移。但事實(shí)上，云的變化方式多種多樣，且目前暫無(wú)可以對(duì)云的變化進(jìn)行精確預(yù)測(cè)的方法，這也是影響DNI預(yù)測(cè)結(jié)果精度的主要原因之一。

4 云層的分類

云層的光學(xué)厚度是決定其對(duì)太陽(yáng)輻射衰減產(chǎn)生影響的主要因素。本文利用XGB算法對(duì)云層的光學(xué)厚度進(jìn)行了分類，具體為：晴空、氣溶膠、薄云、可能為厚云、厚云、藍(lán)色的厚云(下文簡(jiǎn)稱為“藍(lán)厚云”)。由于氣溶膠是影響DNI預(yù)測(cè)結(jié)果精度的主要因素之一，因此將其單獨(dú)作為一類進(jìn)行分類。此外，部分云用肉眼很難判斷是薄云或厚云，因此添加了“可能為厚云” 這一類型。

在提取圖像特征之前對(duì)太陽(yáng)及其鄰近區(qū)域進(jìn)行了插值處理，這樣有助于得到較好的DNI預(yù)測(cè)結(jié)果。

在保證數(shù)據(jù)準(zhǔn)確性的基礎(chǔ)上，本文采用半人工方式提取了一定數(shù)量的云層分類特征樣本，用于云層分類模型的訓(xùn)練，并對(duì)云層分類模型進(jìn)行了驗(yàn)證。不同類型云層的分類結(jié)果占比如表1所示。

表1 不同類型云層的分類結(jié)果占比Table 1 Proportion of classification results of different types of cloud layers

有時(shí)會(huì)存在某些云層分類錯(cuò)誤的情況，比如：可能為厚云被歸類為厚云，藍(lán)厚云被歸類為厚云等，但這是符合預(yù)期的。在這種情況下，云層分類的預(yù)期范圍準(zhǔn)確率可達(dá)90%以上，如表2所示。

表2 云層分類的預(yù)期范圍準(zhǔn)確率Table 2 Accuracy of expected range of cloud classification

將通過(guò)云層分類訓(xùn)練后得到的云層分類模型用于DNI預(yù)測(cè)時(shí)的云層分類，將分類后得到的云層類型作為超短期DNI預(yù)測(cè)模型的特征之一。

5 超短期DNI預(yù)測(cè)

5.1 基于XGB算法的超短期DNI預(yù)測(cè)模型訓(xùn)練

根據(jù)之前計(jì)算得到的云的運(yùn)動(dòng)矢量，利用外推法獲得未來(lái)0～30 min內(nèi)(間隔1 min)云層可能會(huì)遮擋太陽(yáng)的30個(gè)圖像點(diǎn)，然后對(duì)這30個(gè)圖像點(diǎn)進(jìn)行圖像特征提取，并利用云層分類模型完成云層分類。以通過(guò)云層分類提取到的圖像特征和云層類型分類結(jié)果作為超短期DNI預(yù)測(cè)的特征。特征提取過(guò)程完全由軟件自動(dòng)完成。

全天空成像儀旁邊安裝有DNI測(cè)量?jī)x器，用于采集DNI數(shù)據(jù)。受云層和其他因素的影響，DNI會(huì)產(chǎn)生衰減，實(shí)測(cè)DNI的衰減率DDNI可表示為：

式中：DNIt為理論的DNI值；DNIm為實(shí)測(cè)的DNI值。

建立超短期DNI預(yù)測(cè)模型，并將DDNI作為其輸出值。該模型的訓(xùn)練數(shù)據(jù)選用2018～2019年期間收集的所有多云天氣中以天為單位篩選出的云的形態(tài)較為穩(wěn)定時(shí)的數(shù)據(jù)，然后根據(jù)云層類型分類結(jié)果和實(shí)測(cè)DNI的衰減率對(duì)篩選出的數(shù)據(jù)進(jìn)行進(jìn)一步過(guò)濾，過(guò)濾掉云層類型和實(shí)測(cè)DNI的衰減率明顯不符的數(shù)據(jù)，使通過(guò)訓(xùn)練得到的超短期DNI預(yù)測(cè)模型更為穩(wěn)定。最后利用XGB算法完成超短期DNI預(yù)測(cè)模型的訓(xùn)練。

不同云層類型時(shí)DDNI的取值范圍如表3所示。

表3 不同云層類型時(shí)DDNI的取值范圍Table 3 Value range of DDNI for different cloud layers types

5.2 基于XGB算法的超短期DNI預(yù)測(cè)

超短期DNI預(yù)測(cè)時(shí)同樣根據(jù)云的運(yùn)動(dòng)矢量，采用外推法獲得未來(lái)0～30 min內(nèi)(間隔1 min)云層可能會(huì)遮擋太陽(yáng)的30個(gè)圖像點(diǎn)，對(duì)應(yīng)30個(gè)時(shí)刻，并提取這30個(gè)時(shí)刻的特征；然后利用超短期DNI預(yù)測(cè)模型對(duì)這30個(gè)時(shí)刻的DNI進(jìn)行預(yù)測(cè)。

相較于實(shí)測(cè)DNI，利用超短期DNI預(yù)測(cè)模型得到不同天氣(晴天、多云、陰天)下提前10 min時(shí)的預(yù)測(cè)DNI曲線，并將該預(yù)測(cè)曲線與理論DNI曲線和實(shí)測(cè)DNI曲線進(jìn)行對(duì)比，具體如圖12所示。

從圖12中可以看出，不同天氣情況下，超短期DNI預(yù)測(cè)模型得到的提前10 min的預(yù)測(cè)DNI與實(shí)測(cè)DNI曲線較為吻合。

針對(duì)2018～2019年中的152天且每天太陽(yáng)高度角大于5°的時(shí)段，將基于XGB算法的超短期DNI預(yù)測(cè)模型得到的所有1～30 min 的DNI預(yù)測(cè)結(jié)果進(jìn)行統(tǒng)計(jì)分析，具體如表4所示。

圖12 不同天氣下超短期DNI預(yù)測(cè)模型提前10 min的預(yù)測(cè)DNI曲線與理論及實(shí)測(cè)DNI曲線的對(duì)比Fig.12 Comparison between predicted DNI curve of ultrashort-term DNI prediction model 10 minutes ahead of time and theoretical and actual DNI curves under different weather conditions

表4 不同提前預(yù)測(cè)時(shí)下DNI預(yù)測(cè)結(jié)果的MAE與RMSETable 4 MAE and RMSE of DNI prediction results under different advance prediction time

(續(xù)表)

從表4可以看出，不同提前預(yù)測(cè)時(shí)間的總體DNI預(yù)測(cè)值的平均絕對(duì)誤差(MAE)為116.0 W/m2，均方根誤差(RMSE)為203.0 W/m2。在未使用可持續(xù)模型的情況下，從這2個(gè)值來(lái)看，預(yù)測(cè)結(jié)果的精度較高。

將實(shí)測(cè)的DNI數(shù)據(jù)與超短期DNI預(yù)測(cè)模型的預(yù)測(cè)結(jié)果進(jìn)行對(duì)比，當(dāng)以“某個(gè)時(shí)刻的DNI值是否大于或小于400 W/m2”作為判斷標(biāo)準(zhǔn)時(shí)，超短期DNI預(yù)測(cè)模型的預(yù)測(cè)結(jié)果的正確率均大于85%，符合塔式太陽(yáng)能熱發(fā)電站運(yùn)行過(guò)程中對(duì)DNI預(yù)測(cè)的指標(biāo)要求。而通過(guò)進(jìn)一步分析發(fā)現(xiàn)，影響DNI預(yù)測(cè)結(jié)果精度的主要因素為氣溶膠及云的不規(guī)則變化。

6 結(jié)論

本文未使用可持續(xù)模型，而是采用基于XGB算法及全天空成像儀圖像實(shí)現(xiàn)了超短期(0～30 min)的DNI預(yù)測(cè)。通過(guò)使用基于太陽(yáng)位置的全天空成像儀標(biāo)定方法，完成了全天空成像儀的姿態(tài)校正和畸變矯正；準(zhǔn)確擬合了晴空背景參數(shù)并建立晴空庫(kù)，為云層的檢測(cè)提供了基礎(chǔ)；通過(guò)對(duì)太陽(yáng)及其臨近區(qū)域進(jìn)行插值修復(fù)，使1 min內(nèi)的DNI預(yù)測(cè)成為可能。根據(jù)云的運(yùn)動(dòng)矢量，利用外推法得到未來(lái)云層可能會(huì)遮擋太陽(yáng)的圖像點(diǎn)，并提取圖像特征，用于基于XGB算法的超短期DNI預(yù)測(cè)模型的訓(xùn)練，最終以訓(xùn)練后的該模型完成超短期DNI預(yù)測(cè)，并得到了較好的結(jié)果，符合塔式太陽(yáng)能熱發(fā)電站運(yùn)行過(guò)程中對(duì)DNI預(yù)測(cè)的指標(biāo)要求。