樊華, 王文旭, 孫杰, 李曉陽

(河南省地震局, 鄭州 450018)

地震、滑坡、泥石流、山洪等是人類最具威脅的幾種自然災害,尤其是破壞性地震會給國民經濟和生命財產帶來巨大損失。災害發(fā)生后,及時獲取災情信息,如房屋倒塌、道路中斷、橋梁坍塌等狀況,對應急救援、抗震救災以及災后重建等具有重要意義[1]。

隨著遙感和無人機技術的發(fā)展,遙感影像越來越多地應用于救災應急的各個方面。災情發(fā)生后,遙感技術以其獨特的視角提供了災區(qū)的大量圖像信息,如何及時有效處理和應用這些遙感影像是一個值得研究的方向。

利用同一地區(qū)災害前后的遙感影像檢測重要地物的變化,即可獲得大尺度地物(道路、河流、大型水庫等)的變化,也能獲取小尺度地物(如房屋、橋梁等)的變化,為了解災情和救災提供有用的信息。

目前提出的圖像變化檢測方法非常多,從檢測圖像粒度的角度將檢測方法分為3類[2]:基于像素的變化檢測;基于對象的變化檢測;混合變化檢測。

基于像素的變化檢測以像素為處理單位,利用光譜信息來生成差異圖并進行分類。其優(yōu)勢是簡單、快速、有效,應用廣泛;劣勢是對配準精度要求高,否則對應像素之間出現(xiàn)位移偏差,進而產生偽信息。這類方法主要有直接比較法(差分法,比值法等)、變換類方法(主成分分析法,纓帽變化法,變化向量分析法等)、機器學習法[3](支持向量機,決策樹,神經網絡等)。

基于對象的變化檢測通常包含圖像分割、特征計算、地物提取3個方面。其檢測單位不再是像素,而是經圖像分割后得到的有意義的對象(像素簇)。其優(yōu)勢是整合了光譜信息、紋理信息和幾何信息,提高了變化檢測精度;劣勢是現(xiàn)有的分割方法很難針對同一類研究目標實現(xiàn)最佳分割。這類方法主要有直接對象比較法(對象級變化向量分析法,馬爾科夫隨機場,條件隨機場等)、對象類別比較法(對象級分類后比較法,柵格GIS矢量集成法等)以及深度學習法[4]。

混合變化檢測融合了上述兩類方法的特點來進一步降低噪聲和雜散點對變化檢測的影響,進一步提高了檢測精度。這類方法主要有光譜混合分析法、模糊聚類分析法、像素級和對象級混合法等。如朱春宇等[5]提出將深度置信網絡與數(shù)學形態(tài)學融合模型對高分辨率遙感影像進行建筑物變化檢測,有效提高檢測率;靖娟利等[6]基于多特征融合的反向傳播神經網絡對高分影像進行變化檢測;逄濤等[7]基于特征提取網絡,提出一種在時間維度上基于像素位置偏移的圖像特征差異增強方法,在多種數(shù)據(jù)集上均表現(xiàn)出較好的性能。

對精確配準的兩幅遙感圖像做變化檢測,差值法和比值法具有較高的檢測精度和較快的計算速度,而其他方法,如分類法及主成分方法等已被證明效果不如差值法和比值法[8]。鑒于計算速度以及實用程度等因素,現(xiàn)采用基于像素的直接比較法(差值法)并結合基于結構特征的邊緣檢測法對遙感影像的重要地物(建筑物)做變化檢測。

建筑物變化檢測通常多用于城市規(guī)劃與管理、國土資源管理和利用、地理信息以及城市數(shù)據(jù)更新等方面。在震害方面,多是在震后進行建筑物的損毀檢測和災情評估,而對地震前后建筑物的變化檢測工作做的較少,故現(xiàn)提出利用改進的Wv_Canny邊緣檢測方法對地震前后的遙感影像的重要地物(建筑物)進行變化檢測。

現(xiàn)針對時序遙感影像,首先利用三維塊匹配(block matching 3D,BM3D)去除高斯噪聲,之后對待配準圖像利用尺度不變特征變換(scale-invariant feature transform,SIFT)方法進行精確配準,然后求取配準圖像和參考圖像的差分圖像并利用改進的Wv_Canny邊緣檢測方法對其進行邊緣檢測,最后融合配準圖像和參考圖像,并根據(jù)邊緣檢測結果在融合圖像上自動標注災害前后兩圖像的差異區(qū)域。

上述步驟中對差分圖像的邊緣檢測是關鍵一步,使用的檢測方法Wv_Canny在兩方面對經典Canny方法做了改進:一是采用自適應貝葉斯小波閾值去噪和自適應中值濾波替代高斯濾波,使Wv_Canny方法能同時去除高斯噪聲和椒鹽噪聲,而Canny方法僅去除高斯噪聲;二是采用梯度的矢量計算替代梯度的標量計算,使Wv_Canny方法能同時利用亮度信息和色度信息,而Canny方法僅利用亮度信息。改進后的Wv_Canny方法在弱邊緣檢測和邊緣提取的完整連續(xù)性方面較經典Canny方法有了一定的提升,將有利于震害影像中弱邊緣(如淺色建筑物邊緣)的提取,有利于增強災后破壞性地物邊緣的提取效果,更好地服務于地震前后遙感影像的重要地物的變化檢測。

1 遙感影像去除高斯噪聲

噪聲嚴重影響圖像特征點的識別、提取和匹配。故精確配準前需去噪。目前傳統(tǒng)降噪方法中去高斯噪聲效果最好的方法是BM3D[9]。其基本思想是將非局部塊匹配和變換域去噪相結合,算法分兩個階段:基礎估計(閾值去噪)和最終估計(維納濾波),每階段均包含相似塊分組、變換域協(xié)同濾波和重疊塊加權平均3個步驟。

1.1 基礎估計

1.1.1 相似塊分組

(1)

1.1.2 變換域協(xié)同硬閾值濾波

協(xié)同濾波有兩個作用,一是降低相似塊的噪聲,二是增強相似塊的稀疏性,有利于變換域圖像系數(shù)的稀疏表達。過程表達式為

(2)

1.1.3 重疊塊加權平均

由于候選塊Zx會出現(xiàn)在不同參考塊ZxR的3D數(shù)組中,去噪后返到圖像原位置的相似塊會重疊,故需加權平均。對重疊塊加權平均后得到基礎估計階段的去噪結果,表達式為

(3)

1.2 最終估計

1.2.1 相似塊分組

1.2.2 變換域協(xié)同維納濾波

基礎估計的目的除了引導最終估計階段的分組外,還為了求維納濾波器所需的收縮系數(shù),收縮系數(shù)表達式為

(4)

(5)

1.2.3 重疊塊加權平均

與基礎估計類似,最終估計階段的重疊塊加權平均結果為

(6)

圖1給出了BM3D方法去除彩色圖像、灰度圖像和遙感圖像高斯噪聲的效果圖,其中圖1(a)和圖1(b)來自Dabov(2007)。

圖1 BM3D去除彩色圖像、灰度圖像和遙感圖像高斯噪聲結果Fig.1 BM3D removes Gaussian noise from color image, gray image and remote sensing image

2 圖像配準

對同一場景,由于拍攝圖像時的各種環(huán)境條件不同(如地點的移動、光照的明暗、拍攝角度的差異,大氣的變化等),所拍攝圖像會出現(xiàn)不同程度的縮放、旋轉以及灰度值的差異等。圖像配準的主要任務就是尋找最佳的幾何變換關系和灰度變換關系,為消除兩幅或多幅圖像間存在的各種差異進行空間和灰度上的最佳匹配。其數(shù)學表達式為

I2(i,j)=g{I1[f(i,j)]}

(7)

式(7)中:f為二維幾何變換;g為一維灰度變換,反映了成像設備因成像條件不同所對應的灰度變化。i和j分別是圖像寬度和高度的像素坐標。當只有幾何變換時,式(7)則簡化為I2(i,j) =I1[f(i,j)],I1(i,j)和I2(i,j)分別為待配準圖像和參考圖像。

圖像配準是圖像識別、融合、拼接的基礎,在震害目標識別、地物毀壞的檢測、遙感圖像融合等領域有廣泛的應用。通過圖像配準,可以比較不同圖像中的常見特征。借此,可能會發(fā)現(xiàn)河流是如何遷移的,一個區(qū)域是如何被淹沒的,街區(qū)建筑物是如何被震害摧毀的等信息。

配準方法主要分為3類:基于區(qū)域灰度信息的方法、基于變換域信息的方法和基于特征信息的方法[10]。

(1)基于區(qū)域灰度信息的方法是利用灰度計算圖像的相似程度,這類方法根據(jù)圖像特點選擇相似度量函數(shù),在選定的幾何變換模型參數(shù)空間內,按所選搜索算法進行搜索,從而得到相似度最大的幾何變換參數(shù)。優(yōu)點是一般無需對圖像做復雜的預處理,避免由于特征提取而引起的誤差,且容易實現(xiàn);缺點是不能直接用于配準圖像的非線性形變和幾何變形大的圖像,且計算復雜度高。

(2)基于變換域信息的方法有傅氏變換、小波變換、Warsh變換等。常用的是基于傅里葉變換的配準方法。

(3)基于特征的方法是在空間域或變換域中提取那些對縮放、旋轉、平移等具有不變性的特征。在空間域,常用特征有點、線、邊緣、輪廓、區(qū)域中心以及閉合區(qū)域等,如角點、交叉點、曲率間斷點、線段構成的邊緣和輪廓、成片湖泊、森林、建筑群等。這類方法計算復雜度低、穩(wěn)健性好,可適用于部分存在復雜幾何變形的圖像。

目前基于特征點的圖像配準技術,效果好的主要有兩種經典配準算法SIFT和加速穩(wěn)健特征(speeded up robust features,SURF)等[11]。前者有更高的精度,后者有更快的速度,二者都具有良好的穩(wěn)健性[12-13]。其主要區(qū)別是圖像多尺度空間的構建方法:SIFT在金字塔的每一層進行高斯濾波并提取特征點,通過金字塔體現(xiàn)尺度不變性,而SURF使用Hessian矩陣,其特征點的正確匹配數(shù)量多于SIFT,金字塔僅僅用來檢測特征點[14]。

Juan等[15]從耗時和性能(尺度、旋轉、模糊、亮度變化、仿射變換)等6個方面給出了SIFT和SURF的比較,如表1所示。

表1 兩種方法的參數(shù)比較表Table 1 Parameter comparison table of the two methods

從表1中數(shù)據(jù)可知, SURF速度更快,對實時性要求高的工程更適用;SIFT配準精度更高,適合于實時性要求不高而精度要求高的工程;SIFT的性能優(yōu)勢是善于應對圖像的縮放變化和旋轉變化,弱項是不善于應對圖像的亮度變化;SURF的性能優(yōu)勢是善于應對圖像的亮度變化,而不善于應對圖像的旋轉變化。

SIFT算法既具有對旋轉、尺度縮放、亮度變化保持不變性,同時對視角變化、仿射變換、噪聲也保持一定程度的穩(wěn)定性。故使用SIFT算法進行圖像配準,流程圖如圖2所示。

圖2 SIFT算法圖像配準流程Fig.2 SIFT algorithm image registration process

3 Wv_Canny算法邊緣檢測

建筑物、街道等均有其規(guī)則的邊緣特征,但災后(比如震災或洪災后)其規(guī)則的邊緣特征會變得雜亂無章。據(jù)此特征可對災害前后的精確配準后的圖像做圖像差分和邊緣檢測,進而對倒塌或沖毀的建筑物等重要地物進行識別并做出標注,便于災情統(tǒng)計和救助等工作。

針對傳統(tǒng)Canny算法中高斯濾波去噪模糊邊緣的問題,提出了改進方法Wv_Canny(W表示小波閾值去噪和自適應中值濾波,v表示矢量梯度),將多尺度分析、矢量梯度融入經典Canny方法。選擇去噪保邊的小波閾值去噪和中值濾波替代高斯濾波,用矢量場梯度替代Canny方法中的標量場梯度,提高了去噪保邊緣效果,且小波穩(wěn)態(tài)分解的平移不變性保證了小波域中圖像邊緣不發(fā)生偏移。

Wv_Canny算法分三步:首先用小波穩(wěn)態(tài)分解閾值去噪去除高斯噪聲,自適應中值濾波去除椒鹽噪聲;然后計算彩圖的矢量場梯度(灰度圖仍用Canny方法計算梯度);最后進行非極大值抑制、雙閾值邊緣檢測和連接。

3.1 小波穩(wěn)態(tài)分解閾值去噪過程

步驟1選擇B樣條小波和Haar小波為小波基[16],最大分解尺度取3。

步驟2對含噪圖像進行3級小波穩(wěn)態(tài)分解。

步驟3采用自適應貝葉斯閾值作為去噪閾值[17]。

由于VisuShrink閾值對信號“過扼殺”,SureShrink閾值對噪聲“過保留”,而BayesShrink閾值隨空間和尺度變化,可有效減緩“過扼殺”和“過保留”現(xiàn)象,表達式為

(8)

噪聲標準差為

(9)

式(9)中:Yj為含噪圖像j尺度小波分解的對角分量;median表示取中值。

無噪圖像信號標準差為

(10)

其中,噪聲圖像的方差為

(11)

式(11)中:Mj為j尺度高頻子帶的系數(shù)個數(shù)。

步驟4使用漸進半軟閾值函數(shù)對小波系數(shù)去噪

硬閾值去噪在邊緣處易產生振蕩,軟閾值易模糊邊緣。吳安全等[18]提出的漸進半軟閾值函數(shù),可有效避免邊緣處的振蕩和模糊。公式如下。

(12)

步驟5對閾值去噪后的小波系數(shù)重構得到去噪結果。

3.2 自適應中值濾波

中值濾波在椒鹽噪聲密度較高時,通過增大濾波器窗口可改善去噪效果,但會導致細節(jié)模糊。而自適應中值濾波根據(jù)噪聲情況自動調整窗口大小,可兼顧去噪和保細節(jié),故采用自適應中值濾波。

3.3 矢量梯度

研究表明,彩圖邊緣檢測在紅綠藍(red-green-blue,RGB)比色調-飽和度-亮度(hue-saturation-value,HSV)、亮度-色度-濃度 (lightness,u-value, v-value,LUV)等空間能取得更優(yōu)的結果[19],選擇RGB空間。彩圖邊緣檢測算法有兩類:基于單通道的(單通道邊緣合成、單通道梯度合成)和基于矢量的方法[20]。單通道的兩種合成法均沒考慮彩圖分量間的相關性,容易丟失邊緣,而矢量法則可得到更多好的邊緣[21]。下面給出Zenzo[22]的矢量梯度公式。

矢量梯度模的表達式為

(13)

矢量梯度方向的表達式為

(14)

式(14)中:i和j分別為圖像寬度和高度的像素坐標。

gii、gjj、gij分別為RGB空間矢量點積。

(15)

式(15)中:u,v為RGB空間的矢量,表達式分別為

(16)

式(16)中:r、g和b分別為單位矢量。

4 應用實例與結果分析

實際應用通常是給定一幅包含場景和地物的參考圖像,之后在另一幅同場景但地物不同的待配圖像中實現(xiàn)兩者的匹配。兩幅圖像中的場景和地物一般是尺度不同、旋轉角度不同以及亮度不同,這些差異是最常見的情形。故按尺度不變性、旋轉不變性和處理實際遙感圖像的可行性3個方面進行試驗,以驗證本文方法的有效性和處理實際遙感圖像的可行性。

4.1 圖像尺度不變性仿真試驗

為驗證本文方法的有效性,圖3給出了3種不同尺度的配準和檢測結果的對比。圖3有3行6列,每行分別為參考圖像、待配準圖像、配準圖像、差分圖像(參考圖和配準圖的差分)、邊緣檢測及差異標注圖、融合及差異標注圖(參考圖和配準圖融合并標注差異);每列的圖像含義相同,由每列的標題給出。從3行的配準、差分、檢測和標注結果看,3種尺度圖像的結果基本相同,視覺效果幾乎無差別,這說明SIFT配準方法具有良好的尺度不變性。

第一列為三幅相同的參考圖像;第二列為三幅大小不同的待配準圖像;331×354表示圖像的大小;每行6幅圖是一個完整體;每列3幅圖的含義均相同,由每列的標題給出圖3 SIFT算法圖像配準尺度不變性Fig.3 SIFT algorithm image registration scale invariance

4.2 圖像旋轉不變性仿真試驗

為驗證本文方法的有效性,圖4給出了3種不同旋轉角度的圖像配準、差分、邊緣檢測和標注結果的對比。圖4分布結構和含義同圖3。第2列從上到下圖像旋轉角度分別為10°、20°和30°。從3行各圖像的配準、差分、邊緣檢測和標注結果來看,3種旋轉角度的最后結果基本相同,視覺效果幾乎無差別。但從第4列的圖像可看出,隨著旋轉角度的逐漸增大,配準精度逐漸出現(xiàn)了微小的降低,但并不影響參考圖像和配準圖像之間主體差異的檢測和標注,這說明SIFT配準方法具有良好的旋轉不變性。

第一列為三幅相同的參考圖像;第二列為三幅旋轉角和大小均不同的待配準圖像;331×354表示圖像的大小;每行6幅圖是一個完整體;每列3幅圖的含義均相同,由每列的標題給出圖4 SIFT算法圖像配準旋轉不變性Fig.4 SIFT algorithm image registration rotation invariance

4.3 街區(qū)車輛以及簡單地物變化的識別和圈定

為驗證本文方法的可行性和實用性,采用實際無人機圖像進行識別試驗,選取了包含多種地物的影像進行算法測試,圖像來自多旋翼型大疆經緯M200無人機拍攝,如圖5(a)和圖5(b)所示。

兩大綠圈主要圈定了地物和黑色轎車的差異區(qū),兩小綠圈和兩小綠橢圓主要圈定了地物差異區(qū);紅色標注了白色轎車差異處,藍色標注了黑色轎車差異處,黃色標注了地物差異處圖5 街區(qū)圖像中車輛變化的識別和標注Fig.5 Identification and annotation of vehicle changes in block images

圖5為街區(qū)圖像的車輛和簡單地物變化的識別和標注圖。第一行為參考圖像(a)、待配準圖像(b)以及配準圖像(c);第二行為差分圖像(d)、邊緣檢測及標注圖像(e)以及融合及標注圖像(f)。從參考圖像(a)和待配準圖像(b)可知,圖5(b)中道路上少8輛白色轎車和大綠圈內的1輛黑色轎車。從差分圖像(d)中可明顯看到8輛白色轎車的大致形狀,但大綠圈內的黑色轎車由于與路面顏色接近而不易看出,大綠圈內的地物也不易明顯看出。若將大綠圈放大后就可明顯看到黑色轎車和白色轎車旁邊的地物;兩小綠圈處的地物差異可明顯看到;兩個綠色橢圓處的地物差異在放大后亦可明顯看到。圖5(e)和圖5(f)給出了8輛白色轎車(紅色標注)和1輛黑色轎車(藍色標注)以及3處簡單地物差異(黃色標注)的識別和標注。結果表明:本文方法對圖像中同一場景不同時間的兩張圖像中的車輛和簡單地物的變化能夠很好地識別并標注出來,表明了方法的可行性和實用性。

4.4 建筑物變化的識別和范圍圈定

為驗證本文方法對重要地物變化識別的可行性和實用性,采用建筑物較多的實際遙感影像進行識別試驗,影像來自Planet(拍攝時間2022年4月,拍攝地點為上海),含大、中、小型各類建筑物、道路、植被等,如圖6(a)和圖6(b)所示。

圖6 房屋的變化檢測和范圍圈定Fig.6 House change detection and range delineation

圖6為建筑物變化的識別和標注圖。第一行為參考圖像(a)、待配準圖像(b)以及配準圖像(c);第二行為差分圖像(d)、邊緣檢測及標注圖像(e)以及融合及標注圖像(f)。從圖6(a)和圖6(b)可知,有3處差異明顯的地方:①兩圖像左上方紅色框所標區(qū)域,圖6(b)中40多棟房屋變化;②兩圖像右下方藍色框所標區(qū)域,圖6(b)中10多棟房屋變化;③兩圖像右側長條形綠色框所標區(qū)域,沿河的單排建筑物變化。這可從圖6(d)差分圖像中明顯看到這3處的差異區(qū)域,隨后是差分圖像的Wv_Canny雙閾值邊緣檢測結果并對差異區(qū)域做了紅色標注。

通過對比實際建筑物變化面積和本文方法識別標注出的變化面積,正確率79%,誤檢率21%,無漏檢率。結果表明:本文方法對同一場景不同時間兩幅遙感影像中建筑物的變化能夠較好地識別并標注出來,可復制推廣于復雜、重要地物變化的識別,如震害前后建筑物破壞的快速識別等,具有一定的可行性和實用性。

4.5 震后建筑物變化的識別和范圍圈定

影像來自長光衛(wèi)星,圖7(a)拍攝時間為2022年11月6日、圖7(b)拍攝時間為2023年2月7日;拍攝地點為土耳其卡赫拉曼馬拉什省蒂爾克奧盧地區(qū),如圖7(a)和圖7(b)所示。

圖7 地震前后建筑物的變化檢測和范圍圈定Fig.7 Change detection and range delineation of buildings before and after earthquake

圖7為建筑物變化的識別和標注圖。圖7(a)為參考圖像(地震前拍攝)、圖7(b)為配準圖像(地震后拍攝)、圖7(c)為參考圖像和配準圖像融合后的差異標記圖像。通過對比圖7(a)和圖7(b),在圖7(b)中分別用紅線畫出了震后建筑物未倒塌的3個范圍,用綠線畫出了震后建筑物倒塌的兩個范圍,用藍線畫出了震前震后地物有變化的兩個范圍。圖7(c)是用本文方法經過配準(SIFT圖像配準方法)、差分、邊緣檢測(Wv_Canny雙閾值邊緣檢測方法)、差異標注等步驟得到的結果。

通過對比震前震后建筑物的倒塌面積和非地震因素造成的變化面積(如震前人工搭建或拆除一些簡單房屋以及其他地物等造成的變化)以及本文方法識別標注出的倒塌面積和變化面積,其正確率為78%,誤檢率為22%。應用實例結果表明:本文方法具有一定的可行性和實用性,可應用于震害前后建筑物等破壞的快速識別。

5 結論

噪聲會嚴重影響圖像特征點的識別、提取和匹配,BM3D是目前傳統(tǒng)降噪方法中去高斯噪聲效果最好的方法,不僅有高的信噪比,還能有效保護圖像的邊緣信息,視覺效果也很好。這將有利于震害影像中高斯噪聲的消除和圖像邊緣的保護,為震害影像的后續(xù)處理如提取滑坡泥石流邊界、建筑物邊緣等工作奠定基礎。

SIFT和SURF圖像配準算法是目前基于特征點效果好的兩種圖像配準技術,前者精度高,適合于精度要求高而實時性要求低的工程,應對圖像的縮放變化和旋轉變化是其強項;后者計算速度快,適合于實時性要求高的工程,應對圖像的亮度變化是其強項;二者都具有很好的穩(wěn)健性。

Wv_Canny邊緣檢測方法將多尺度分析、矢量梯度融入經典Canny方法,能檢測出震害影像中更多弱邊緣,如滑坡泥石流邊界、淺色建筑物邊緣等,其保邊作用提升了邊緣提取的完整連續(xù)性,使邊緣檢測結果整體效果得到了提升。

仿真實驗表明,SIFT配準方法具有良好的尺度不變性和旋轉不變性;應用實例表明,通過對遙感影像的精確配準、圖像差分、邊緣檢測、融合并標注變化前后圖像的差異區(qū)域,實驗結果按建筑物變化面積比較,正確率為78%～79%,誤檢率為21%～22%,無漏檢率。仿真和實際實驗結果表明了本文方法對于遙感影像中簡單地物和復雜地物變化的識別均具有可行性和實用性,對災后倒塌和毀壞的建筑物等重要地物差異區(qū)域能有效識別并標注,有利于災情的及時了解和救援等工作。