李世偉,王召巴,楊建生
(1.中北大學(xué) 信息與通信工程學(xué)院,山西 太原030051;2.美國(guó)鮑爾州立大學(xué) 地理系,曼西 印第安那47036;3.中北大學(xué) 山西省光電信息與儀器工程技術(shù)研究中心,山西 太原030051)
近年來(lái),我國(guó)城市進(jìn)程加快,地物分布信息變化頻繁,大量人工構(gòu)建面成為了城市的主體,造成諸如城市熱島效應(yīng)、城市內(nèi)澇等問(wèn)題時(shí)有發(fā)生。作為人類活動(dòng)頻繁的城市區(qū)域,其典型地物主要包括植被、水體、不透水面(包括道路、房屋等)和裸土區(qū),對(duì)這些地物分布信息的及時(shí)掌握,是有效監(jiān)測(cè)城市生態(tài)環(huán)境和制定規(guī)劃決策的關(guān)鍵。隨著科學(xué)技術(shù)的不斷發(fā)展,基于遙感的城市地物信息獲取方法得被廣泛應(yīng)用[1-4]。
目前,被廣泛認(rèn)可并能得到較好提取結(jié)果的方法是利用地物光譜特征進(jìn)行組合運(yùn)算的指數(shù)法。每一個(gè)指數(shù)都可以針對(duì)一類地物信息進(jìn)行專題提取,通過(guò)波段運(yùn)算和閾值劃分,可以快速獲得專題信息。這些指數(shù)中,比較知名的包括歸一化植被指數(shù)(NDVI)[5]及其改進(jìn)型指數(shù)[6]、歸一化水體指數(shù)(NDWI)和改進(jìn)型歸一化水體指數(shù)(MNDWI)[7]、歸一化不透水層指數(shù)(NDISI)[8,9]。但值得注意的是,這些指數(shù)運(yùn)算結(jié)果只能對(duì)一種專題信息進(jìn)行顯示,無(wú)法使用戶直接獲得多類地物的分類信息,因此,本文就這一問(wèn)題展開(kāi)研究。
為了信息收集和對(duì)實(shí)驗(yàn)結(jié)果進(jìn)行驗(yàn)證的方便,本文選取了作者所在城市太原市。太原主城區(qū)位于太原河谷平原北部,屬盆地結(jié)構(gòu),地處太行山與呂梁山之間,屬暖溫帶大陸性季風(fēng)氣候,季節(jié)性明顯。冬季通常寒冷干燥,到了春季氣溫迅速回升,雨量集中在夏季,晝夜溫差較大;秋季天高氣爽,氣溫變化起伏較小,全年光熱資源充足,年均溫度介于7.6-10.2℃之間,年均降水量為420-460毫米,年均降水量滿足作物的基本生長(zhǎng)需求,但降水量的季節(jié)差異較大,易造成干旱災(zāi)害和水土流失。主城區(qū)自北向南,汾河水從城區(qū)中部流過(guò),在城區(qū)汾河中段的西面,分布著城區(qū)內(nèi)唯一的湖泊-晉陽(yáng)湖。
近年來(lái),太原市經(jīng)歷了城市擴(kuò)張后再綠化的過(guò)程,土地覆蓋類型變化劇烈。為了對(duì)預(yù)想的方法進(jìn)行驗(yàn)證,選取了一景LANDSAT-TM 影像,該影像的成像時(shí)間為2011年4月19 日(軌道:PATH-125,ROW-33),影像中心經(jīng)度為112.311度,緯度36.046度。在這個(gè)月份,研究區(qū)的耕地作物還沒(méi)有長(zhǎng)出地面,因此對(duì)應(yīng)區(qū)域全部為土壤覆蓋,而除去農(nóng)作物外的其它植被諸如樹(shù)木和草地等已經(jīng)開(kāi)始生長(zhǎng),呈現(xiàn)出了較為明顯的植被特征。此外,為了對(duì)最后的分類結(jié)果精度進(jìn)行驗(yàn)證,選取了相同時(shí)相的一景SPOT 高分辨率影像作為參考。
為了利用指數(shù)運(yùn)算的快速性優(yōu)勢(shì),并同時(shí)獲得多類地物信息,按照以下步驟進(jìn)行實(shí)驗(yàn):首先,使用指數(shù)法對(duì)專題信息進(jìn)行提??;而后,針對(duì)每個(gè)指數(shù)的提取圖像結(jié)果進(jìn)行閾值劃分和二值化處理;最后,在不破壞專題信息特征的基礎(chǔ)上,利用偽彩色賦色的方法直接對(duì)專題信息的二值化圖層進(jìn)行圖層堆疊和偽彩色賦色,以顏色作為分類依據(jù),直接呈現(xiàn)多類地物分布信息。即:①使用指數(shù)分別對(duì)水體、植被和不透水層專題信息進(jìn)行提??;②利用高分辨率影像,通過(guò)目視解譯和對(duì)比對(duì)信息模糊地點(diǎn)進(jìn)行辨別,得到每個(gè)專題信息圖層的合理閾值,并使用二值化方法得到專題信息;③對(duì)得到的圖層進(jìn)行偽彩色合成實(shí)驗(yàn),得到可以反映四類地物信息的最好色彩組合方式;④以高分辨率影像為參考,得到僅通過(guò)顏色信息獲得的地物分類精度結(jié)果,以驗(yàn)證該方法的可行性。
所選TM 影像和SPOT 影像均為二級(jí)產(chǎn)品,已經(jīng)過(guò)系統(tǒng)級(jí)輻射校正和幾何校正。為了獲得更為精確的結(jié)果,首先對(duì)影像進(jìn)行了幾何精校正,在采用多項(xiàng)式模型和雙線性插值法進(jìn)行了幾何校正和圖像重采樣后,從而得到了精校正影像。其次,使用了1:25萬(wàn)山西省行政區(qū)邊界矢量數(shù)據(jù)生成了太原市主城區(qū)的矢量邊界,利用此邊界生成研究區(qū)AOI數(shù)據(jù)層對(duì)影像進(jìn)行裁剪。最終獲得的太原市主城區(qū)影像如圖1右圖所示。
圖1 太原市主城區(qū)遙感影像
由于遙感圖像在成像時(shí),會(huì)受到成像角度、太陽(yáng)高度角等因素的影響,所獲取的圖像數(shù)據(jù)均以灰度值進(jìn)行顯示,會(huì)造成 “同物異譜”和 “同譜異物”現(xiàn)象的發(fā)生,為了在一定程度上消除這些影響,最為有效的辦法就是將灰度圖像反演為反射率影像。為此,首先將研究區(qū)的TM 影像的第1-5波段和第7波段轉(zhuǎn)換為反射率影像。其過(guò)程如下
(1)將6個(gè)波段的圖像DN 值轉(zhuǎn)換為輻射亮度值
式中:R——像元的輻射亮度值,DN 原始影像中的像元數(shù)字灰度,Gain、Bias——USGSLANDSAT-TM 用戶手冊(cè)中提供的對(duì)應(yīng)波段的增益數(shù)據(jù)和偏置數(shù)據(jù)(見(jiàn)表1)。
表1 LANDSAT-TM 增益和偏置數(shù)據(jù)
(2)反射率計(jì)算
式中:ρ——像元的反射率,D——日地天文單位距離(通常取1),R——式(1)中的像元輻射亮度值,ESUNI——大氣頂層太陽(yáng)輻射平均值(見(jiàn)表2),——太陽(yáng)天頂角。計(jì)算天頂角的計(jì)算公式為
式中:SUN_ELEVATION——對(duì)應(yīng)影像頭文件中的太陽(yáng)高度角,在所使用的影像頭文件中,其值為57°。
表2 大氣頂層太陽(yáng)輻射平均值
根據(jù)以上的計(jì)算過(guò)程,利用ERDAS IMAGINE 遙感圖像處理軟件的MODLE MAKER 建立了相應(yīng)的計(jì)算模型并進(jìn)行計(jì)算,得到了反射率反演結(jié)果影像,結(jié)果如圖2所示(RGB=432)。
為了在將來(lái)使用本方法時(shí)能夠提高工作效率,利用ERDAS IMAGINE遙感圖像處理軟件的MODLE MAKER建立了相應(yīng)的指數(shù)計(jì)算模型,從而在對(duì)其它區(qū)域或本區(qū)域其它時(shí)相使用時(shí),快速生成對(duì)應(yīng)的專題信息圖層。
通過(guò)歸一化植被指數(shù)NDVI進(jìn)行植被信息的獲取,其表達(dá)式為
式中:ρNIR——第4 波段的反射率,ρRed——第三波段的反射率,其結(jié)果如圖3所示。其中灰度值相對(duì)較高的區(qū)域?yàn)橹脖粎^(qū)域,灰度值較低的為非植被區(qū)域。為了使結(jié)果更具有可讀性,在結(jié)合高分辨率影像進(jìn)行了多次植被信息提取的閾值選取后,對(duì)NDVI進(jìn)行了二值化處理,處理閾值為0.374,大于閾值區(qū)域顯示為白色,反之為黑色,二值化處理圖像如圖4所示。
圖3 NDVI處理結(jié)果
水體信息的提取采用了由徐涵秋教授于2005年提出的改進(jìn)型歸一化水體指數(shù)(MNDWI),其表達(dá)式為
式中:Green——綠光波段,MIR——中紅外波段,分別為L(zhǎng)ANDSAT-TM 的第2 波段和第5 波段,提取結(jié)果如圖5所示,其中灰度值較高處為水體。同樣為了使結(jié)果更具有可讀性,在結(jié)合高分辨率影像進(jìn)行了多次植被信息提取的閾值選取后,對(duì)MNDWI進(jìn)行了二值化處理,處理閾值為0.197,大于閾值區(qū)域顯示為白色,反之為黑色,二值化處理圖像如圖6所示。
圖4 NDVI處理結(jié)果
圖5 MNDWI閾值處理結(jié)果
圖6 MNDWI閾值處理結(jié)果
不透水層信息的提取,目前已經(jīng)有多種方法。這些方法中,大部分都是基于Ridd于1995 年提出的VIS 模型。但是這些模型通常不能體現(xiàn)水體信息,因此,在不透水層提取前需要先對(duì)水體信息進(jìn)行掩膜操作,將水體信息提出后在進(jìn)行下一步的操作,這種做法勢(shì)必會(huì)將誤差引入不透水層提取的工作中來(lái)[10]。2008年徐涵秋教授在總結(jié)和比較已有指數(shù)的基礎(chǔ)上,重新通過(guò)對(duì)地物各自在影像上的光譜特征以及地物之間的光譜特征差異分析,找到了能夠使地物差異最大化的波段進(jìn)行組合計(jì)算,提出了歸一化不透水層指數(shù)(NDISI)來(lái)對(duì)不透水層信息進(jìn)行提取。其表達(dá)式如下
式中:TIR——熱紅外波段,對(duì)TM 數(shù)據(jù)就是第六波段。NIR、MIR1——TM第4和第5波段。本文首先使用了NDISI指數(shù)對(duì)影像進(jìn)行處理,但經(jīng)過(guò)對(duì)提取結(jié)果進(jìn)行解譯和實(shí)地驗(yàn)證,我們?cè)贜DISI提取結(jié)果中發(fā)現(xiàn):太原武宿機(jī)場(chǎng)等新建不透水面的灰度值高于水體,而大量道路、樓宇等建筑的灰度出現(xiàn)了低于水體灰度的情況,甚至出現(xiàn)少量不透水層與水體灰度值相等的情況。
因此,即使不考慮少量不透水層與水體灰度值相等的情況,單純使用NDISI提取不透水層也需要兩個(gè)以上的閾值才可以將不透水層信息提取出來(lái),并且除了每增加一個(gè)閾值會(huì)相應(yīng)增加誤差外,后續(xù)還需要對(duì)處理結(jié)果進(jìn)行掩膜、合并等大量操作?;谶@種考慮,我們采取了另外的方法進(jìn)行了不透水面的提取??紤]到指數(shù)方法基于地物光譜特征,我們重新對(duì)所使用影像中的典型地物進(jìn)行光譜特征采集(如圖7所示),圖中縱坐標(biāo)為T(mén)M 影像波段編號(hào),縱坐標(biāo)為多個(gè)典型地物樣點(diǎn)的灰度均值。通過(guò)對(duì)該特征曲線組進(jìn)行分析后,得到了以下結(jié)果:
圖7 四類地物光譜特征
不透水層在第3波段和第4波段的差異均大于其它3種地物包括水體、植被和土壤,若將4波段和3波段相除,可達(dá)到增大地物灰度差異的目的;不透水層在5 波段和6波段的光譜特征與其它3種地物呈現(xiàn)出相反的情況,即不透水層第5波段灰度大于其第6波段,而其它3種地物變現(xiàn)為第6波段灰度大于其地物波段;在第4波段和第5波段,水體也呈現(xiàn)出于其它3種地物的相反特征。
基于以上結(jié)果,并考慮到第6波段與其它6個(gè)波段分辨率不一致,且MNDWI的提取結(jié)果中水體呈現(xiàn)了較好的提取效果,因此采取了以下方案對(duì)不透水面進(jìn)行提?。?/p>
(1)首先對(duì)影像的第4和第3波段進(jìn)行比值操作(RI,式(7)),對(duì)原圖像的不透水層和水體信息同時(shí)進(jìn)行增強(qiáng),結(jié)果如圖8所示。
通過(guò)式(7)計(jì)算后,不透水層信息和水體信息同時(shí)得到增強(qiáng),且水體灰度略高于高于不透水層,并且所用不透水層的灰度值大小相當(dāng),沒(méi)有出現(xiàn)NDISI中出現(xiàn)的不透水層灰度既有大于水體的部分又有小于水體灰度的情況的發(fā)生。
圖8 RI結(jié)果影像
(2)由于在RI圖層中,水體和不透水層均具有較高的灰度,在MNDWI中的水體信息得到增強(qiáng)而不透水層信息被抑制,為此,首先將RI圖層在進(jìn)行閾值選取后的二值化處理,處理閾值為1.159,而后,將該二值化圖層與MNDWI二值化圖層進(jìn)行差值運(yùn)算后,最終得到了不透水層的信息圖層(如圖9所示)。
為了實(shí)現(xiàn)分類,對(duì)所獲取的專題圖層進(jìn)行分析,考慮到對(duì)NDVI、MNDWI、RISI減MNDWI的3個(gè)提取結(jié)果已經(jīng)進(jìn)行了二值化處理,又考慮到ERDAS軟件可以快速實(shí)現(xiàn)對(duì)3個(gè)圖層進(jìn)行堆疊和偽彩色賦色的操作,因此直接對(duì)3個(gè)二值化圖層進(jìn)行了圖層堆疊,并分別賦予不同的偽彩色,這樣做的原因是:
在植被專題信息、水體專題信息和不透水層專題信息層的二值化圖像中,其各自的專題信息區(qū)域灰度值為1,而其它信息區(qū)域灰度為0,因此,若將3個(gè)專題信息圖層賦予不同的偽彩色并合成后,各專題信息在合成結(jié)果中的顏色應(yīng)具有較大差異。而且,由于土壤在3個(gè)專題圖層上的灰度值都為零,無(wú)論對(duì)哪一圖層賦予除黑色外的偽彩色后,土壤都將呈現(xiàn)為黑色。因此,在使用堆疊三圖層并偽彩色賦色后,四類地物應(yīng)該可以從顏色進(jìn)行區(qū)分。
圖9 不透水層提取結(jié)果
基于以上考慮,在將3 個(gè)圖層進(jìn)行圖層堆疊操作后,對(duì)NDVI二值化圖層賦為綠色,不透水面二值化圖層賦為紅色,水體二值化圖層賦為藍(lán)色(如圖10所示)。
為了對(duì)本文得到的結(jié)果進(jìn)行驗(yàn)證,得到該方法的可行程度,在最后得到的偽彩色合成結(jié)果中隨機(jī)選取了200個(gè)點(diǎn),在結(jié)合高分辨率影像進(jìn)行對(duì)比后,獲得的具體分類精度見(jiàn)表3。
表3 分類精度評(píng)價(jià)
可以看出,該方法可得到83%的總體分類精度,為了使評(píng)價(jià)結(jié)果更具有客觀性,通過(guò)式 (8)對(duì)分類結(jié)果進(jìn)行了KAPPA 系數(shù)的計(jì)算。
式中:N——參與精度統(tǒng)計(jì)的像元數(shù)量,Xii——分類正確的像元數(shù)量Xi+是表格中第i行像元的數(shù)量,X+i——表格中第i列像元的數(shù)量,所得KAPPA 系數(shù)結(jié)果為0.803。
圖10 偽彩色合成結(jié)果
通過(guò)參考同時(shí)期的高分辨率影像,本文使用的方法所造成的誤分情況多發(fā)生在汾河南端的河道和沿岸發(fā)生,通過(guò)實(shí)地考察和詢問(wèn),并結(jié)合相關(guān)資料分析發(fā)現(xiàn),在影像獲取時(shí)段內(nèi),汾河太原段中斷和南段水流較小,河道裸露,小型死水面較多,其面積小于影像中一個(gè)像元所對(duì)應(yīng)的面積,造成了遙感影像中混合像元的出現(xiàn),因此,在一定程度上影響了分類精度。
采用建模方式建立了多個(gè)指數(shù)模型,并通過(guò)運(yùn)算得到了相應(yīng)的專題信息,采用閾值劃分和二值化處理的方式后,在沒(méi)有使用任何分類算法的情況下,利用偽彩色合成方法同時(shí)得到了包括水體、植被、不透水面和裸土共四類地物的分布信息。其中,利用NDVI和MNDWI分別對(duì)植被和水體進(jìn)行了提取,結(jié)合不透水層和水體兩者在MNDWI圖層以及經(jīng)過(guò)光譜分析計(jì)算建立的RI圖層中的特點(diǎn),對(duì)二值化后的RI圖層和MNDWI圖層進(jìn)行差值運(yùn)算后得到了不透水層信息。在將得到的分類結(jié)果圖像與高分辨率參考圖像對(duì)比后,得到了分類結(jié)果的總體分類精度為83%,Kappa系數(shù)0.803。結(jié)果表明,該方法通過(guò)色彩區(qū)分就可直接快速得到研究區(qū)四類地物的分類結(jié)果和分布影像,省去了常規(guī)分類方法中必需的類別合并、分類賦色等過(guò)程,模型實(shí)現(xiàn)容易,閾值改變僅需在模型中改變數(shù)值即可,實(shí)現(xiàn)了四類城市典型地物分布信息的快速獲取,可為地物分布信息的快速監(jiān)測(cè)提供參考技術(shù)手段。
[1]ZHAO Ming.Analysis on the evolvement of Wenzhou City’s morphological elements [J].Urban Studies,2009,16 (7):28-32 (in Chinese).[趙明.溫州城市形態(tài)構(gòu)成要素演化研究[J].城市發(fā)展研究,2009,16 (7):28-32.]
[2]SONG Yang,WANG Shijun,YANG Yanru,et al.Study on the urban internal spatial structure in northeast China [J].Areal Research and Development,2009,28 (4):18-23 (in Chinese).[宋飏,王士君,楊艷茹,等.東北地區(qū)城市內(nèi)部空間結(jié)構(gòu)研究 [J].地域研究與開(kāi)發(fā),2009,28 (4):18-23.]
[3]MU Fengyun,ZHANG Zengxiang,TAN Wenbin.Analysis on the spatial-temporal characteristics of Guangzhou City’s spatial morphologic evolution [J].Geo-Information Science,2007,9 (5):94-98 (in Chinese). [牟鳳云,張?jiān)鱿?,譚文彬.廣州城市空間形態(tài)特征與時(shí)空演化分析 [J].地球信息科學(xué),2007,9 (5):94-98.]
[4]LIU Yaxuan,ZHANG Xiaolei,LEI Jun,et al.Spatial expansion and driving forces of oasis cities in Xinjiang,China [J].Journal of Desert Research,2011,31 (4):1015-1021 (in Chinese).[劉雅軒,張小雷,雷軍,等.新疆綠洲城市空間擴(kuò)展特征及其驅(qū)動(dòng)力分析 [J].中國(guó)沙漠,2011,31 (4):1015-1021.]
[5]SONG Dongmei,ZHANG Qian,YANG Xiuchun,et al.Spatial and temporal characteristics of MODIS vegetation index in the source region of three rivers on Qinghai-Tibet plateau in China[J].Geographical Research,2011,30 (11):2067-2075(in Chinese). [宋冬梅,張茜,楊秀春,等.三江源區(qū)MODIS植被指數(shù)時(shí)空分布特征 [J].地理研究,2011,30 (11):2067-2075.]
[6]ZHU Gaolong,LIU Yibo,JU Weimin,et al.Evaluation of topographic effects on four commonly used vegetation indices 4[J].Journal of Remote Sensing,2013,17 (1):210-234 (in Chinese).[朱高龍,柳藝博,居為民,等.4種常用植被指數(shù)的地形效應(yīng)評(píng)估 [J].遙感學(xué)報(bào),2013,17 (1):210-234.]
[7]XU Hanqiu.Comment on the enhanced water index(EW I):A D iscussion on the creation of awater index [J].Geo-Information Science,2008,10 (6):776-780 (in Chinese).[徐涵秋.從增強(qiáng)型水體指數(shù)分析遙感水體指數(shù)的創(chuàng)建 [J].地球信息科學(xué),2008,10 (6):776-780.]
[8]XU Hanqiu.A new remote sensing index for fastly extracting impervious surface information [J].Geomatics and Information Science of Wuhan University,2008,33 (11):1150-1153 (in Chinese).[徐涵秋.一種快速提取不透水面的新型遙感指數(shù)[J].武漢大學(xué)學(xué)報(bào)·信息科學(xué)版,2008,33 (11):1150-1153.]
[9]XU Hanqiu.Analysis of impervious surface and its impact on urban heat environment using the normalized difference impervious surface index(NDISI) [J].Photogrammetric Engineering &Remote Sensing,2010,76 (5):557-565.
[10]MU Fengyun,ZHANG Zengxiang,TAN Wenbin.Analysis on the spatial-temporal characteristics of Guangzhou City’s spatial morphologic evolution [J].Geo-Information Science,2007,9 (5):94-98 (in Chinese).[牟鳳云,張?jiān)鱿?,譚文彬.廣州城市空間形態(tài)特征與時(shí)空演化分析 [J].地球信息科學(xué),2007,9 (5):94-98.]