胡學(xué)謙,孫 林,隋淞蔓,顧曉鶴,孫 宇,何玉龍
(1.山東科技大學(xué) 測繪與空間信息學(xué)院,山東 青島 266590;2.國家農(nóng)業(yè)信息化工程技術(shù)研究中心,北京 100097)
云層的存在,阻擋了太陽輻射能量的傳輸,嚴(yán)重影響了遙感圖像的利用率,給地物分類、地表參數(shù)定量反演等應(yīng)用帶來較大的困難,但在遙感成像中,云是難以避免的,為了減少云的影響,需要對云進(jìn)行檢測。由于雪和云在可見光和近紅外波段具有相似的反射率特征,當(dāng)遙感影像中存在積雪時,易導(dǎo)致云檢測精度下降[1]。
目前主要的云檢測方法有閾值法、人工神經(jīng)網(wǎng)絡(luò)法和統(tǒng)計學(xué)方法。閾值法因為具有算法簡單、速度快和檢測精度高的優(yōu)點而被廣泛使用[2]。閾值法云檢測的基本思想是基于云和晴空像元的反射率差異和亮度溫度差異進(jìn)行波譜分析和特征通道選擇,使用單通道、多通道組合的反射率、亮溫等設(shè)置一系列云檢測閾值。常用的閾值法有ISCCP方法[3]、APOLLO法[4]、CLAVR[5]法和CO2薄片法[6]等。傳統(tǒng)閾值法使用典型地物和云的光譜分析或者基于圖像統(tǒng)計得出閾值,沒有充分考慮下墊面和云形態(tài)的復(fù)雜性,云檢測精度較低。國內(nèi)外對積雪識別也進(jìn)行了廣泛研究。歸一化差值積雪指數(shù)(NDSI)的方法由于其簡單且精度較高而被廣泛使用[7],并在此基礎(chǔ)上發(fā)展出了增強(qiáng)型積雪指數(shù)(ENDSI)法[8]。1995年,Hall等人在NDSI方法的基礎(chǔ)上提出了SNOWMAP算法[9],用于使用MODIS數(shù)據(jù)繪制全球積雪地圖,是目前最常用的積雪識別方法之一。以上這些方法都沒有與云檢測有效結(jié)合起來。
先驗數(shù)據(jù)庫支持的云檢測閾值自動生成(CDAG)算法是一種針對不同衛(wèi)星傳感器的云檢測閾值生成方法[10]。算法通過在高光譜分辨率的AVIRIS數(shù)據(jù)上人工標(biāo)注的云像元和晴空像元,模擬不同傳感器的云和晴空像元庫,自動、充分分析云和地表不同組分的光譜特征差異,找到最合適的云檢測閾值,具有檢測效率高、適應(yīng)性強(qiáng)和檢測效果好的優(yōu)點。由于CDAG算法是基于AVIRIS高光譜數(shù)據(jù)建立的像元庫運行的,受飛行平臺的限制,像元庫中雪像元較少,并且在構(gòu)建像元庫時,忽略了云和雪反射率特征的相似性,導(dǎo)致對存在雪覆蓋的影像云檢測精度下降。
為改進(jìn)原算法在云檢測中對雪像元識別效果差,易誤判的問題,本文提出了一種增強(qiáng)雪識別的改進(jìn)CDAG算法。采用Landsat8陸地成像儀(Operational Land Imager,OLI)數(shù)據(jù)進(jìn)行驗證,結(jié)合多種云檢測算法,利用改進(jìn)后的像元庫進(jìn)行迭代計算,自動生成云檢測閾值,實現(xiàn)增強(qiáng)雪識別的云檢測。經(jīng)驗證明,算法能夠準(zhǔn)確識別出影像中的云、雪像元。
為了驗證增強(qiáng)雪識別的改進(jìn)CDAG算法的有效性,需要研究區(qū)域同時存在積雪和云層覆蓋,而青藏高原地區(qū)平均海拔4 000 m以上,部分地區(qū)終年積雪,地貌類型簡單,以巖石、凍土、冰川和裸地為主,并且太陽輻射強(qiáng),大氣條件純凈,適合進(jìn)行增強(qiáng)云雪識別云檢測算法的驗證,是本研究的良好研究區(qū)域[11]。選擇青藏高原地區(qū)78°E~85°E,29°N~36°N區(qū)域作為研究區(qū),所選區(qū)域地處青藏高原西部,覆蓋喀喇昆侖山、帕米爾高原、昆侖山西麓和喜馬拉雅山西麓,年積雪覆蓋日數(shù)(SCD)60~120天,平均雪深10 cm左右,冬春季節(jié)(11月—次年5月)云量45%~70%,盛行高云[12]。研究區(qū)內(nèi)存在不同類型的云層和積雪場景,存在大量混合像元,增加了云雪識別的難度,是檢驗算法的理想?yún)^(qū)域。
CDAG云檢測算法的基礎(chǔ)在于高光譜像元庫的構(gòu)建,像元庫是利用AVIRIS數(shù)據(jù)構(gòu)建的。如表1所示,AVIRIS是第一個測量370~2 507 nm太陽輻射光譜的成像光譜儀,具有光譜分辨率、信噪比和空間分辨率高的特點,并通過224個連續(xù)的光譜通道以10 nm為間隔測量上行輻射[13],能夠探測到地物在光譜特性上更微小的差異[14]。
表1 AVIRIS傳感器特征
本研究選擇Landsat8 OLI數(shù)據(jù)作為實驗數(shù)據(jù)。Landsat8衛(wèi)星于2013年2月在美國加州發(fā)射,其搭載的OLI共有11個波段,刈幅寬度為185 km,重訪周期為16 d,全色波段的空間分辨率為15 m,其余波段為30 m[15]。云和雪的反射波譜特征差異在波長0.38~2.5 μm的可見光波段和近紅外波段最為明顯,也是云檢測的常用波段[16-17],實驗數(shù)據(jù)選擇的Landsat8 OLI 1~7波段和第9波段即位于這一波段范圍內(nèi)。表2為選擇的波段及其波長范圍。
表2 選用的Landsat8 OLI波段及波長范圍
選擇研究區(qū)域冬春季節(jié)的6景Landsat8 OLI數(shù)據(jù)進(jìn)行驗證,所選數(shù)據(jù)均存在不同類型積雪及云覆蓋。數(shù)據(jù)如表3所示。
表3 驗證數(shù)據(jù)信息
本文選用的數(shù)據(jù)為L1T級數(shù)據(jù)產(chǎn)品,已經(jīng)過預(yù)處理,但為了保證數(shù)據(jù)的一致性,進(jìn)行實驗驗證前,對數(shù)據(jù)做統(tǒng)一的輻射定標(biāo)[18]。使用影像元數(shù)據(jù)中的輻射校正系數(shù),將DN值轉(zhuǎn)化為大氣上層反射率。
定標(biāo)公式為:
(1)
式中,ρλ為大氣上層表觀反射率;Mρ,Aρ為從Landsat8 OLI元數(shù)據(jù)中獲得的定標(biāo)系數(shù);Qcal為Landsat8 OLI原始影像DN值;θsz為太陽天頂角。
CDAG算法是一種自動閾值的云檢測算法。云檢測的基本原理是根據(jù)云在可見光、近紅外和短波紅外波段與其他下墊面(植被、土壤、水體、巖石和建筑等)的光譜特征差異對云進(jìn)行檢測[19]。CDAG算法使用預(yù)先標(biāo)識出云和典型地物的高光譜像元庫,利用數(shù)據(jù)模擬的方式來實現(xiàn)不同衛(wèi)星的云檢測閾值確定,其實現(xiàn)主要包括以下幾個部分:
(1) 高光譜像元庫的構(gòu)建。在選擇AVIRIS高光譜數(shù)據(jù)建立云像元庫和晴空像元庫時,為了讓云檢測算法更加完善,構(gòu)建的閾值有廣泛適應(yīng)性,選擇的云和地表類型的種類要全面。云像元庫要同時包含厚云、薄云和碎云等像元;晴空像元庫要同時包括裸地、巖石、建筑、水體和林地等像元,并且要有足夠的像元數(shù)量。
(2) 利用高光譜像元庫模擬不同的多光譜傳感器的像元。CDAG算法的模擬方法采用了Lian等人[20]提出的方法。根據(jù)波譜響應(yīng)函數(shù)將窄波段的AVIRIS像元庫數(shù)據(jù)進(jìn)行模擬得到寬波段的待測數(shù)據(jù)像元庫,模擬公式如下[21]:
(2)
(3) 基于模擬得到的多光譜傳感器的云和晴空像元庫,進(jìn)行統(tǒng)計分析,采用機(jī)器學(xué)習(xí),自動計算不同波段和波段組合在不同閾值下的云像元識別正確率和晴空像元誤判率。閾值確定原理如圖1所示。
圖1 閾值確定原理Fig.1 Principle of threshold determination
為了確定每種算法的最佳閾值,將云像元庫和晴空像元庫中的像元數(shù)據(jù)作為訓(xùn)練樣本,統(tǒng)計反射率的最大值和最小值,使閾值在2個最值之間以步長0.01變化,云像元正確率和地表誤判率也隨之變化,當(dāng)云像元的識別正確率達(dá)到最高,且晴空像元誤判率開始隨閾值變化明顯增加時,對應(yīng)閾值為此種算法的最佳閾值。
CDAG算法在達(dá)到較高云檢測效率和精度的同時,也存在沒有充分考慮云和雪反射率特征的相似性,容易對雪像元產(chǎn)生誤判的問題,導(dǎo)致影像存在積雪時云檢測精度降低。
針對CDAG算法存在的問題,本文提出了增強(qiáng)雪識別的改進(jìn)CDAG算法。云與雪的識別主要是利用云雪在可見光和近紅外波段存在的反射率差異。表4為利用AVIRIS高光譜數(shù)據(jù)繪制的400~2 500 nm波段范圍的云、雪、巖石和裸地等研究區(qū)常見目標(biāo)的光譜曲線。
表4 研究區(qū)典型目標(biāo)像元
由表4可以看出,云在950,1 150 nm處有2個明顯的吸收谷,在1 700,2 200 nm處有2個明顯的反射峰,與雪在對應(yīng)波段的反射率存在明顯的差異,可以利用這些波段對云和雪進(jìn)行識別,而巖石和裸土與前二者反射率差異明顯,不易產(chǎn)生混淆算法的改進(jìn)包括以下幾部分:
(1) 利用AVIRIS數(shù)據(jù)構(gòu)建的像元庫,模擬待檢測的Landsat8 OLI數(shù)據(jù)的云像元庫和晴空像元庫[22-24]。與高光譜的AVIRIS傳感器相比,Landsat8 OLI傳感器光譜分辨率較低且每一個波段都覆蓋了AVIRIS的多個波段。通過數(shù)據(jù)模擬將對應(yīng)的若干個連續(xù)的窄波段高光譜數(shù)據(jù)合成對應(yīng)的Landsat8 OLI傳感器的寬波段多光譜數(shù)據(jù),得到Landsat8 OLI數(shù)據(jù)的云像元庫和晴空像元庫。
(2) Landsat8 OLI像元庫的優(yōu)化。由于AVIRIS數(shù)據(jù)構(gòu)建的高光譜像元庫并未充分考慮雪對云檢測的影響,為提高云檢測中云雪識別的效果,選擇研究區(qū)內(nèi)經(jīng)過預(yù)處理的Landsat8 OLI數(shù)據(jù),在影像上選取總數(shù)大于10萬個的不同地表覆蓋的積雪、碎雪等具有代表性的雪像元,提取其在不同波段下的表觀反射率,并將獲取的反射率數(shù)據(jù)加入模擬得到的Landsat8 OLI晴空像元庫,在不同波段組合的云檢測算法下以機(jī)器學(xué)習(xí)的方式自動生成改進(jìn)后的最佳云檢測閾值。
(3) 基于Landsat8 OLI云像元庫和晴空像元庫的反射率差異,對其進(jìn)行波段運算、統(tǒng)計分析、自動生成閾值以及重合度判斷等處理,自動生成云檢測算法。本文主要選取了單波段、雙波段組合、波段比值和波段差值作為云檢測算法的公式,最終的云檢測算法如式(3)~式(6)所示:
ρi>Ti,
(3)
ρi>Ti&ρj>Tj,
(4)
Ta<ρi/ρj (5) Tc<ρi-ρj (6) 式中,i,j為1~7,9;ρi,j代表第i,j波段的反射率;Ti,j代表第i,j波段的最佳閾值;Ta,b代表不同波段比值算法對應(yīng)的最佳閾值;Tc,d代表不同波段差值算法對應(yīng)的最佳閾值。不同波段組合下生成的最佳閾值如表5所示。 表5 不同算法的最佳閾值 (4) 云檢測結(jié)果的生成。利用生成的算法對待測數(shù)據(jù)進(jìn)行實驗,每一種算法都會生成相應(yīng)的云檢測結(jié)果,利用式(7)對不同算法得到的云檢測結(jié)果加權(quán)合成。云概率計算公式為: (7) 式中,F(xiàn)i為云檢測結(jié)果;Qi為對應(yīng)云檢測結(jié)果的權(quán)重;N為云檢測算法總數(shù)。通過加權(quán)合成,得到最終的云檢測結(jié)果,并對檢測結(jié)果進(jìn)行精度驗證和分析。 為了評價增強(qiáng)雪識別的改進(jìn)CDAG算法對于Landsat8 OLI數(shù)據(jù)的云檢測效果,本文對實驗數(shù)據(jù)的云檢測精度進(jìn)行定量分析。采用目視解譯的方式,在測試數(shù)據(jù)的標(biāo)準(zhǔn)假彩色影像上區(qū)分出云覆蓋區(qū)域和雪覆蓋區(qū)域并將其矢量化,與云檢測結(jié)果進(jìn)行比較,分別統(tǒng)計識別正確和識別錯誤的云像元數(shù)量和晴空像元數(shù)量,將像元分為表6中的4類。 表6 像元分類 采用以下幾個指標(biāo)對云檢測精度進(jìn)行評價:云像元識別正確率(CRA)、云像元漏判率(CRM)、晴空像元識別正確率(SRA)和晴空像元漏判率(SRM),4個指標(biāo)的計算如下: (8) (9) (10) (11) 在以上4個指標(biāo)中,CRA和SRA可以反映增強(qiáng)雪識別的改進(jìn)CDAG算法的可靠性;CRM指漏判的云像元在總的云像元中的占比,反映了算法低估云的程度;SRM指漏判的晴空像元在總晴空像元中的占比,反映了算法高估云的程度。 圖2為增強(qiáng)雪識別的改進(jìn)CDAG算法的云檢測結(jié)果示例影像,其中左圖為Landsat8 OLI數(shù)據(jù)5,4,3波段合成的假彩色影像,右圖為對應(yīng)的云雪識別結(jié)果,白色為云,灰色為雪,黑色為下墊面。通過對比左圖的假彩色影像和右圖的檢測結(jié)果,可以看出原數(shù)據(jù)中的云和雪都被幾乎完整地提取出來,并賦予不同顏色加以區(qū)分,體現(xiàn)出算法具有較好的識別效果。 圖2 云雪識別結(jié)果示意Fig.2 Result of cloud and snow identification 根據(jù)云的幾何形態(tài)和密度,可以將云劃分為厚云、薄云和碎云[25]。當(dāng)云層呈大面積的塊狀,分布均勻,且無法透過下墊面信息,則認(rèn)為該云層為厚云。厚云主要包括低云和中云中較厚的云層,一般呈亮白色或灰白色,形狀大小沒有規(guī)律。薄云與厚云的區(qū)別在于能否從影像中透過云層觀察到下墊面信息。由于薄云輕薄和半透明的特性,其光譜信息往往與地表混合。碎云指面積較小,厚度不均勻,在影像中呈現(xiàn)離散狀分布的小云塊,由于體積較小,在云檢測時常常會被遺漏。不同類型云的反射率特征存在差異,本研究對不同云層類型覆蓋條件對改進(jìn)CDAG算法的適應(yīng)性進(jìn)行分析。 圖3為改進(jìn)CDAG算法對不同積雪覆蓋區(qū)域上空厚云的檢測結(jié)果,左圖為5,4,3波段合成的標(biāo)準(zhǔn)假彩色影像;中圖為對云和雪的識別效果圖,白色為云,灰色為雪,黑色為除雪像元外其他晴空像元;右圖為云檢測結(jié)果圖,每種云層條件分別選取4個不同區(qū)域和尺度進(jìn)行分析。 (a) 小尺度區(qū)域 厚云是最常見最典型的云層,與其他地物反射率差異大,與積雪的區(qū)別也相對明顯,較容易識別,識別正確率較高。圖3(a)、圖3(b)為較小尺度下存在厚云的區(qū)域,可以看到,積雪上空的云層可以被完整地識別,積雪也可以較為準(zhǔn)確地識別,云和雪具有清晰的邊界,極少出現(xiàn)誤判。圖3(c)、圖3(d)為較大尺度區(qū)域,可以看到,檢測結(jié)果中云被準(zhǔn)確地識別,圖3(d)中云層空隙漏出的積雪地表也被準(zhǔn)確識別,并在云檢測結(jié)果中剔除,體現(xiàn)出方法的有效性。 圖4為改進(jìn)CDAG算法對Landsat8 OLI數(shù)據(jù)不同積雪覆蓋區(qū)域上空薄云的檢測結(jié)果。薄云具有可以部分顯示下墊面的特點,在常用云檢測波段的反射率較厚云更低,在存在積雪的區(qū)域更容易產(chǎn)生誤判。 (a) 小尺度區(qū)域 圖4(a)、圖4(b)為存在薄云的較小尺度積雪覆蓋區(qū)域,可以看出,云和雪被準(zhǔn)確地識別,極少出現(xiàn)誤判。圖4(c)、圖4(d)為存在積雪和其他地物類型的大尺度區(qū)域,從檢測結(jié)果可以看到,絕大部分積雪和薄云均被準(zhǔn)確識別,但也存在少量小塊積雪被誤判為云。在圖像中薄云和積雪地表的過渡部分,存在一定程度的誤判,這是由于薄云可以部分透過下墊面信息的特點,無法與下墊面劃定明確的邊界,對最終的云檢測精度產(chǎn)生負(fù)面影響。 碎云在影像上一般表現(xiàn)為小面積圖斑,相互獨立且形狀不規(guī)則,同時具有厚云和薄云的特點,反射率變化較大且邊緣較復(fù)雜。圖5為碎云覆蓋區(qū)域的云檢測結(jié)果。 (a) 小尺度區(qū)域 圖5(a)、圖5(b)為存在碎云且有大面積積雪覆蓋的較小尺度影像,可以從中看出,碎云和雪均被較為準(zhǔn)確地提取,較少出現(xiàn)誤判。圖5(c)、圖5(d)為存在碎云和積雪且尺度相對較大的影像,在這2組圖上可以較為清晰地看到碎云和雪的邊界,碎云被較好地識別,受到雪的影響較小,誤判較少。圖5(b)中碎云相對圖5(c)較薄,對比2幅圖像可以直觀看出,圖5(b)的云檢測結(jié)果中出現(xiàn)較多被誤判為云的白色亮點,而圖5(c)檢測結(jié)果更為準(zhǔn)確,這是由于當(dāng)碎云較薄時,體現(xiàn)出薄云的特點,更易產(chǎn)生混合像元,增加了誤判,降低了云檢測的精度,說明對碎云的檢測更為困難。 為了定量評價改進(jìn)云檢測算法的準(zhǔn)確性,對其進(jìn)行精度驗證。計算不同云雪覆蓋條件下云檢測結(jié)果的4個指標(biāo),結(jié)果如表7所示。 表7 云檢測結(jié)果定量評價 表7是圖3、圖4和圖5對應(yīng)的存在積雪覆蓋的厚云、薄云和碎云3類影像云檢測結(jié)果的定量評價和總體的精度評定。從表7可以看到,改進(jìn)的云檢測算法對于Landsat8 OLI數(shù)據(jù)的云檢測達(dá)到了較高的精度,總體的云像元識別正確率和晴空像元識別正確率均達(dá)到了90%以上。參考圖3、圖4和圖5的云檢測結(jié)果影像可以看到,改進(jìn)的算法較少出現(xiàn)云像元和雪像元的誤判,并且對其他晴空像元的誤判率也較低。分別分析厚云、薄云和碎云的定量評價結(jié)果可以看到,不同類型的云對檢測結(jié)果存在影響,厚云的整體檢測精度較高且檢測結(jié)果較穩(wěn)定,云識別的平均正確率達(dá)到了95.2%,晴空像元的平均漏判率為4.6%;而薄云和碎云的檢測精度相對低,云識別的平均正確率分別為91.9%和92.2%,晴空像元的平均漏判率分別為3.8%和4.6%??傮w的云識別平均正確率為93.1%,總體的晴空像元平均漏判率為4.3%。通過精度驗證可以看出,本文提出的增強(qiáng)雪識別的改進(jìn)云檢測算法對Landsat8 OLI數(shù)據(jù)的云檢測結(jié)果具有很高的精度,并且能夠正確區(qū)分云雪像元。 云的存在阻礙了太陽輻射能量的傳輸,增加了地物識別和地表參數(shù)反演的難度,導(dǎo)致遙感圖像利用率下降,為了降低云的影響,需要進(jìn)行云檢測[26]。而在高海拔地區(qū)和中高緯度地區(qū)的寒冷季節(jié),地表常會存在積雪,雪和云的反射波譜特征具有相似性,增加了此類區(qū)域的云檢測難度,對遙感影像的判讀和使用產(chǎn)生消極影響[27]。為了提高對存在積雪覆蓋區(qū)域的云檢測精度,本文提出了一種增強(qiáng)雪識別的改進(jìn)CDAG算法,并對這種方法的有效性進(jìn)行了驗證。本文利用青藏高原地區(qū)的多景Landsat8 OLI影像數(shù)據(jù)對改進(jìn)CDAG算法進(jìn)行了實驗,并對厚云、薄云和碎云3種典型云覆蓋條件下的積雪覆蓋區(qū)域分別進(jìn)行了云檢測。 從圖3、圖4和圖5的云檢測結(jié)果可以看出,本文方法針對不同云層類型均能得到準(zhǔn)確的云檢測結(jié)果,能夠正確地識別出云和雪,在云檢測結(jié)果中云和雪及其他地貌類型具有清晰的邊界,對于云的識別完整、清晰自然,達(dá)到了很好的效果。從表7的精度評價結(jié)果可見,總體的云像元識別正確率達(dá)到了93.11%,晴空像元的識別正確率達(dá)到了95.65%,較傳統(tǒng)的云檢測方法有了很大的提升,顯示了本文方法的有效性。但不同類型云覆蓋下的檢測精度仍存在差異,其中厚云的云像元識別正確率達(dá)到了95.2%,明顯高于薄云和碎云的92%左右,對其進(jìn)行分析,有以下原因:① 厚云的形態(tài)較為穩(wěn)定,云層較厚且與地面有明確的邊界,不易產(chǎn)生混合像元;② 薄云能夠部分反映地面信息,碎云形態(tài)破碎邊界雜亂,這些特性導(dǎo)致了混合像元的產(chǎn)生,從圖中也可以看到這些特征,混合像元的存在模糊了云與雪的波譜特征差異,導(dǎo)致薄云和碎云的檢測精度較厚云偏低。 本文在CDAG算法的基礎(chǔ)上改進(jìn)得到了增強(qiáng)雪識別的改進(jìn)CDAG算法,主要改進(jìn)包括以下2個方面:① 針對雪覆蓋地表,擴(kuò)充了像元庫,增加雪像元在晴空像元中的占比,通過機(jī)器學(xué)習(xí)提高了云檢測閾值的準(zhǔn)確性,提高了對雪像元的識別正確率;② 采用了單波段、多波段組合、波段比值和波段差值多種算法,提高了云檢測的精度。改進(jìn)后的算法在保留了CDAG算法識別正確率高、自動生成閾值等優(yōu)點的基礎(chǔ)上,降低了對雪像元的誤判率,在實驗中取得了良好的檢測效果,在積雪覆蓋的區(qū)域平均云像元識別正確率達(dá)到93.11%。目前算法存在的在薄云和碎云條件下,云檢測精度較厚云條件下降以及像元庫數(shù)據(jù)由人工選擇組成,主觀性較強(qiáng),可能導(dǎo)致數(shù)據(jù)庫的構(gòu)建存在偏差的問題,是下一步需要研究的內(nèi)容。2.3 精度驗證
3 結(jié)果與分析
3.1 厚云的檢測
3.2 薄云的檢測
3.3 碎云的檢測
3.4 精度評定
4 討論
5 結(jié)束語