陈启浩1
,
聂宇靓1,2
,
李林林1
,
刘修国1
|
收稿日期: 2016-11-10; 优先数字出版日期: 2017-11-01
基金项目: 国家自然科学基金(编号：41301477，41471355)；中国地质调查局地质调查项目(编号：DD20160060)
第一作者简介: 陈启浩(1982— )，男，讲师，研究方向为极化SAR信息提取。E-mail：cugcqh@163.com
通讯作者简介: 聂宇靓(1993— )，女，博士研究生，研究方向为极化SAR灾害评估。E-mail：yuliangnie@sina.com
中图分类号: TP753
文献标识码: A
|
摘要
准确地获知灾区的建筑物损毁程度能为抗震救灾和灾后重建提供决策依据。利用震后极化合成孔径雷达(SAR)数据，该文提出了一种综合利用极化分解后多纹理特征的震后建筑物损毁评估方法。首先，用Pauli分解的π/4偶次散射分量剔除非建筑区；其次，用Pauli分解的π/4偶次散射分量的方差特征、对比度特征和Pauli分解的奇次散射分量的对比度特征识别倒塌建筑物，并分别基于区块计算建筑物损毁指数；最后，综合3个纹理特征完成建筑物的损毁评估。采用玉树震后RADARSAT-2数据和东日本大地震后ALOS-1数据的实验验证了所提方法对建筑物损毁评估的有效性，该方法对玉树城区和日本石卷城区的重度、中度和轻度损毁建筑评估的总体精度分别为74.39%和80.26%。与其他方法的对比实验表明，该方法能减少取向角的影响，对存留有少数与方位向平行的完好建筑物的倒塌区、大取向角的完好建筑区的评估更为准确。
关键词
极化SAR, 震害评估, 建筑物, Pauli分解, 灰度共生矩阵
CHEN Qihao1
,
NIE Yuliang1,2
,
LI Linlin1
,
LIU Xiuguo1
Abstract
Accurately obtaining building damage situation in disaster areas provides decision-making support for disaster relief and post-disaster reconstruction. The existing evaluation methods using polarization featuresare greatly influenced by the orientation angle of buildings, and the intact buildings with large orientation anglesare easily misclassified for the collapsed buildings. The combination of multiple features is an effective way to improve the accuracy of building damage assessment, in which the texture features play an important role. This paper proposed a new building damage assessment method which applies the texture features of polarization decomposition components. There are three key procedures in this method. Firstly,the non-building areasare removed by thefiltered π/4 double-bounce scattering component of Pauli decomposition. Secondly,thestrategy of extracting collapsed buildings is built.For building areas, the Variance and Contrast features of Gray-Level Co-occurrence Matrix (GLCM) are calculated on the Pauli decomposition parameters. Specifically, the Variance and Contrast textures are obtained by the π/4 double-bounce scattering component of Pauli decomposition; the Variance texture is given by the odd scattering component of Pauli decomposition. Then, the appropriatethresholds for distinguishing collapsed buildings from intact buildings are determined by these texture features.When the texture feature value of the pixel is smaller than the threshold, then the pixel is classified as collapsed building, otherwise, the pixel is identified as intact building. At last,building damage assessment is implemented. The building damage indexes computed by the three texture features are averaged to get the final building damage index, and then the result of building damage assessment is obtained. This method was validated on RADARSAT-2 fine-mode polarimetric SAR imagery from the Yushu earthquakeacquired on April 21, 2010. The ground-truth map of the buildings in Yushu city, interpreted on high-resolution optical data, was used as a reference for comparison. Compared with the H-α-ρ method, the method with homogeneous texture of span image and the method with integrated Normalized Circular-pol Correlation Coefficient (NCCC) and homogeneous texture feature, the building damage assessment accuracy of the proposed method was the highest, and its overall accuracy was improved to 74.39%.Particularly, the detection rate for serious damage buildings is 100%, and the false alarm rate is 15.62%. In addition, the intact buildings with large orientation angles in northeastern corner of the city were also extracted correctly. Based on the experimental analysis, the main factor of affecting the accuracy is the similarity of Pauli polarization features, the Variance and Contrast texture features betweenthe flat and dense intact building areas and serious damaged building areas. In addition, the effectiveness of the proposed method was verified by the ALOS-1 data after the 2011east Japan earthquake and tsunami, and its overall accuracy was improved to 80.26%. A damage assessment method based on texture features from GLCM of Pauli decomposition components is proposed in this paper. In this method, Pauli decomposition parameter was utilized to remove non-building areas, which was capable of removing the non-building areaseffectively, such as rivers, roads, bare ground and so on. By the comprehensiveutilization of the GLCM texture features of three Pauli decomposition components, the collapsed buildings can be extracted more accurately. By comparison with other methods, the results confirm the validity of the proposed method.
Key words
polarimetric SAR, earthquake damage assessment, building, pauli decomposition, gray-level co-occurrence matrix
1 引 言
地震是破坏性最强的自然灾害之一，经常造成严重的人员伤亡和经济损失。快速、准确地获知震后建筑物损毁情况能及时为抗震救灾提供决策依据。遥感技术具有快速、综合、宏观等优点，在获知震区损毁情况方面具有很大优势(范一大 等，2016)。合成孔径雷达(SAR)以其全天时、全天候、穿透性强、大区域等对地观测特点(郝洪美 等，2012)，成为获取震害信息的一种实用选择。特别是地震发生后经常伴随恶劣的天气条件，因此研究基于SAR数据的地震灾区建筑物损毁评估具有十分重要的意义。
近年来，国内外学者在SAR数据震区建筑物损毁评估方面进行了研究。Yonezawa和Takeuchi(2001)指出地震前后SAR影像的强度相关性和干涉相干特征能一定程度地反映建筑物的损毁程度，并提出了归一化相干系数差的倒塌建筑物探测方法。张景发等人(2002)通过分析张北地震前后倒塌建筑物的相关性、平均灰度差异性和灰度方差差异性，定量提取了建筑物损毁信息。Matsuoka和Yamazaki(2004)以日本神户地震为例，利用震前震后强度影像的相关系数和后向散射系数差构造了一个新的判别指数，并提取倒塌建筑物。Gamba等人(2006)提出了基于边缘特征变化检测和基于像素变化检测的震害信息提取方法。Dekker(2011)用COSMO-SkyMed和TerraSAR-X数据对海地地震进行分析，利用SAR数据归一化强度差和相关系数分别进行了基于网格的平均震害检测和单栋建筑的详细震害检测。Chen和Sato(2013)利用ALOS极化SAR数据分析了2011年东日本大地震和海啸的损毁情况，提出了两个确定损毁等级的指标：震后与震前模型分解的偶次散射分量之比，震后与震前取向角差值的标准差。
上述方法取得了较好的应用效果，但是对数据要求较高，需要同时获取震前震后多时相SAR数据。由于合适的震前数据不易获得，针对震后单时相SAR数据的建筑物损毁评估方法被提出。Guo等人(2009)针对汶川震后的极化SAR数据利用圆极化相关系数ρRRLL、偶次散射分量Pd和各向异性度A提取倒塌建筑物。郭华东等人(2010)针对玉树震区RADARSAT-2极化数据，提出利用H和α参数去除非建筑物，利用圆极化相关系数ρRRLL区分倒塌与完好建筑物的H-α-ρ方法。Zhang等人(2015)用最优极化对比度增强法检测玉树震后RADARSAT-2极化数据中的损毁建筑物，并证明对有一致取向角的建筑区的效果较好。Zhai等人(2016)、Zhai和Huang(2016)分别引入归一化的去取向前后偶次散射分量差、去取向前后偶次散射分量和体散射分量相对贡献变化率之差等参数，修正Wishart监督分类结果，以减少取向角建筑和损毁建筑的误分，从而提高建筑物损毁评估的精度。基于震后SAR数据的建筑物损毁评估方法主要利用极化特征，其评估结果为地震灾后应急和重建提供了重要的科学依据和决策支持。但也存在一些不足：评估结果受建筑物取向角的影响较大，容易将大取向角的完好建筑物误分为倒塌建筑物；对存留有少数与方位向平行的完好建筑物的倒塌区，也易评估不准确。
为解决这些问题，震区建筑物损毁评估方法不再以单一的极化特征为主，而是向多类特征相结合的方向发展。Zhao等人(2013)用H-α-Wishart分类去除非建筑区，综合规范化圆极化相关系数(NCCC)和总功率灰度共生矩阵(GLCM)的同质性纹理特征，评估玉树地震后建筑物的损毁程度，并用高分辨率机载极化SAR数据验证了其有效性。Shi等人(2015)利用极化、干涉和纹理等181个特征，用随机森林分类器对震后机载SAR数据进行分类和建筑物损毁评估，指出纹理特征对倒塌建筑物提取及损毁评估最有效。由此可见，多类特征结合是提高震后建筑物损毁评估精度的一个有效办法，且纹理特征在其中能发挥重要作用。
因此，本文将目标分解的极化特征与灰度共生矩阵纹理相结合，提出一种利用极化分解后多纹理特征的震区建筑物损毁评估方法。首先利用Pauli分解的极化特征剔除非建筑区；然后利用基于Pauli分解极化特征计算的3种GLCM纹理特征分别提取倒塌建筑物，最后综合3种提取结果完成建筑物损毁评估；并利用玉树震后RADARSAT-2极化数据和东日本大地震后ALOS-1极化数据进行验证。
2 本文方法描述
建筑物倒塌前后其散射机制和纹理均会发生变化。考虑到建筑物密集的城区更多表现为低熵散射，选用Pauli相干分解提取极化特征，并结合经典的GLCM纹理描述方法，综合利用震区目标散射机制和纹理特征的差异以准确评估震后建筑物损毁程度。利用极化分解后多纹理特征的建筑物损毁评估方法的流程图如图1所示。
该方法的主要步骤如下：首先，对震区极化SAR数据进行Pauli分解，基于滤波去噪后的π/4偶次散射功率ω′剔除非建筑区。其次，基于Pauli分解的π/4偶次散射功率ω提取灰度共生矩阵的方差(variance)纹理和对比度(contrast)纹理，基于Pauli分解的奇次散射功率u提取灰度共生矩阵的contrast纹理。然后，利用这3种纹理特征单独提取倒塌建筑物并计算区域建筑损毁指数。最后，将基于3种纹理特征计算的区域损毁指数取平均值得到最终损毁评估结果。其中，α，β，γ分别表示基于3个纹理特征区分完好建筑物和倒塌建筑物的阈值。
2.1 Pauli分解与非建筑区剔除
Pauli分解将散射矩阵S分解为多个Pauli基矩阵的复数形式的加权和，每个Pauli基矩阵对应一种基本散射机制(Lee和Pottier，2009)，表示形式如下
|
$\begin{aligned}{S} = & \left[ {\begin{array}{*{20}{c}}{{S_{\rm{HH}}}} & {{S_{\rm{HV}}}}\\[3pt]{{S_{\rm{VH}}}} & {{S_{\rm{VV}}}}\end{array}} \right] = \frac{a}{{\sqrt 2 }}\left[ {\begin{array}{*{20}{c}}1 & 0\\[3pt]0 & 1\end{array}} \right] + \frac{b}{{\sqrt 2 }}\left[ {\begin{array}{*{20}{c}}1 & 0\\[3pt]0 & { - 1}\end{array}} \right] +\\ & \frac{c}{{\sqrt 2 }}\left[ {\begin{array}{*{20}{c}}0 & 1\\[3pt]1 & 0\end{array}} \right] + \frac{d}{{\sqrt 2 }}\left[ {\begin{array}{*{20}{c}}0 & { - j}\\[3pt]j & 0\end{array}} \right]\end{aligned}$
|
(1) |
式中，SHV表示发射电磁波为垂直极化，接收电磁波为水平极化；a，b，c，d均为复数，其值由下式给出
|
$\begin{array}{l}a = \displaystyle\frac{{{S_{{\rm{HH}}}} + {S_{{\rm{VV}}}}}}{{\sqrt 2 }}\;\;\;b = \displaystyle\frac{{{S_{{\rm{HH}}}} - {S_{{\rm{VV}}}}}}{{\sqrt 2 }}\\c = \displaystyle\frac{{{S_{{\rm{HV}}}} + {S_{{\rm{VH}}}}}}{{\sqrt 2 }}\;\;\;d = j\displaystyle\frac{{{S_{{\rm{HV}}}} - {S_{{\rm{VH}}}}}}{{\sqrt 2 }}\end{array}$
|
(2) |
Pauli分解将目标的散射过程相干分解为4种散射机制：平坦表面的奇次散射，方向角为0°的角反射器产生的二面角散射(偶次散射)，方向角为45°的角反射器产生的二面角散射(π/4偶次散射)，以及散射矩阵的所有不对称分量。
在满足互易定理的单站情况下，SHV=SVH，使得d=0，此时Pauli分解可简化为3个基矩阵。奇次散射、偶次散射和π/4偶次散射分量的功率依次可描述为
|
$\begin{array}{l}u = 10{\rm{lg}}\left( {{{\left| a \right|}^2}} \right)\\[7pt]v = 10{\rm{lg}}\left( {{{\left| b \right|}^2}} \right)\\[7pt]\omega = 10{\rm{lg}}\left( {{{\left| c \right|}^2}} \right)\end{array}$
|
(3) |
一般地，非建筑区主要包括河流、道路及裸地等目标，这些地物以奇次散射为主，其他类型散射均较弱；而完好建筑物的屋顶及倒塌建筑物的奇次散射也较强。建筑物倒塌后偶次散射降低，奇次散射增强；但完好建筑和倒塌建筑的π/4偶次散射均明显强于非建筑物区。
图2为玉树震后的RADARSAT-2极化SAR图像中完好建筑物、倒塌建筑物和非建筑物区样本Pauli分解后各分量的均值。可以看出，非建筑区和建筑区(包括完好建筑物和倒塌建筑物)在π/4偶次散射分量的差异最大。因此，根据震后灾区地物的散射机制，本文利用Pauli分解的π/4偶次散射分量功率ω来区分非建筑区和建筑区。
为了降低相干斑噪声的影响，本文利用3×3窗口均值滤波后的π/4偶次散射分量功率ω′剔除非建筑物。通过统计分析确定区分建筑区与非建筑区的阈值ωt，将ω′<ωt的像素分为非建筑区，ω′≥ωt的像素分为建筑区。
阈值的确定方法如下：首先选择两类样本统计其特征值的分布，以样本特征值的重叠区间为阈值所处范围[t1, t2](t1表示特征值较大类别的样本点最小值，t2表示特征值较小类别的样本点最大值)；然后在该阈值范围内搜索，取样本点区分精度最高时对应的阈值为最终阈值。以RADARSAT-2数据划分建筑区和非建筑区的阈值ωt确定为例，建筑区和非建筑区的样本在特征ω′上的特征值分布如图3所示。通过统计分析确定阈值所处范围为[–14.85, –13.15]，分别如图3中的橙色和蓝色直线所示。当阈值为–14.85时，建筑区样本的检出率为100%，有部分非建筑区的样本被识别为建筑区；当阈值为–13.15时，非建筑区样本的检出率为100%，有部分的建筑区样本被识别为非建筑区。在[–14.85, –13.15]范围内搜索，选取对建筑区和非建筑区的区分精度最高的值–13.5作为最终的分割阈值，如图3中的黑色虚线所示。
2.2 利用Pauli极化分量GLCM纹理的倒塌建筑提取
GLCM由Haralick等人(1973)提出的一种通过图像灰度的空间相关性来描述纹理的方法。通过统计图像中一定距离和一定方向上的两个像元灰度之间的相关性，描述图像的纹理信息。其定义如下：在图像范围内，统计所有从灰度级为i的像元点，沿θ方向移动距离D到达灰度为j的像元点的概率p(i, j)，形成的方阵即GLCM。其中距离D可分解为水平方向的Dx和垂直方向的Dy，θ通常取0°，45°，90°，135° 4个方向。
GLCM反映了图像灰度关于方向、相邻间隔及变化幅度的综合信息。利用GLCM可分析图像的局部模式和排列规则，为了更直观地描述纹理，一般用灰度共生矩阵的二阶统计量来表示。本文主要使用方差和对比度两种二阶统计量表示的纹理特征。
(1) 方差(variance)纹理特征
|
$\begin{array}{l}{\rm{variance}} = \mathop \sum \limits_{i = 1}^N \mathop \sum \limits_{j = 1}^N {\left( {i - \mu } \right)^2}p\left( {i,j} \right)\\\mu = \mathop \sum \limits_{i = 1}^N \mathop \sum \limits_{j = 1}^N i \cdot p\left( {i,j} \right)\end{array}$
|
(4) |
式中，N为图像灰度级数，μ表示GLCM的均值。方差纹理特征反映了像元值与其周围局部区域的均值的偏差。当图像中灰度变化较大时，GLCM方差特征值较大。
(2) 对比度(contrast)纹理特征
|
${\rm{contrast}} = \mathop \sum \limits_{i = 1}^N \mathop \sum \limits_{j = 1}^N {\left( {i - j} \right)^2}p\left( {i,j} \right)$
|
(5) |
对比度能反映灰度局部变化的情况，灰度差大的像素越多，对比度值越大，表明灰度局部变化剧烈。
地震发生后，重度损毁区由于建筑物倒塌严重，形成大片废墟，在雷达影像上其相邻像元所包含地物几乎都为废墟，像元间差异较小，整体较为均匀，导致重度损毁区的方差和对比度纹理特征均表现出低值。完好建筑区及中轻度损毁区地物分布较为复杂，包括高度不一的建筑，建筑和废墟的混杂等，相邻像元包含地物差异较大，使得这些区域的方差和对比度纹理特征均表现出高值。
图4为RADARSAT-2数据中完好建筑区和倒塌建筑区样本(如图5(a)所示)在π/4偶次散射功率ω的variance纹理、contrast纹理和奇次散射功率u的contrast纹理上的特征值统计图。可以看出，总体上完好建筑物的方差和对比度纹理特征值都大于倒塌建筑物的特征值。因此本文采用阈值法，基于Pauli分解极化特征计算的3种GLCM纹理特征提取倒塌建筑物。
2.3 区域建筑震损评估
基于提取的倒塌建筑物结果，本文用建筑物损毁指数来描述该区域内建筑物的损毁程度。损毁指数越大，建筑损毁程度越严重。区域的建筑物损毁指数定义为
|
$Cr = {N_{\rm{c}}}/{N_{\rm{b}}}$
|
(6) |
式中，Nc表示该区域范围内倒塌建筑物像素数目，Nb表示该区域内所有建筑物像素数目。
为了充分利用多种极化参数的纹理特征，本文将3种纹理特征分别计算的建筑物损毁指数平均后得到综合的区域倒塌率。最后根据各区域的综合倒塌率划分对应的损毁程度。参照现有的划分方法(Zhai 等，2016；Zhao 等，2013)，确定建筑物损毁程度的划分规则为：Cr≤30%为轻度损毁，30%<Cr≤50%为中度损毁，Cr>50%为重度损毁。
3 实验与分析
本文中的实验，包括纹理特征的提取、建筑损毁程度的评估以及精度评价，都是通过IDL8.3编程实现；纹理特征的样本统计分析利用Matlab2012实现。程序运行环境为64位win7操作系统。
3.1 RADARSAT-2数据实验
为了验证本文方法的有效性，实验采用2010年玉树地震后的RADARSAT-2精细模式极化数据，并与H-α-ρ方法、总功率同质性纹理的方法、综合NCCC与灰度共生矩阵同质性纹理的方法进行定性和定量的对比。
2010年4月14日7时49分，青海省玉树藏族自治州玉树县发生里氏7.1级大地震，城区大量房屋倒塌造成了极大的生命和财产损失。本文所用实验数据为玉树城区2010年4月21日右视升轨的RADARSAT-2精细模式极化数据，影像大小为560×270像素，空间分辨率约8 m，入射角为21°，方位向进行3视平均去噪处理。图5(a)为实验数据的PauliRGB合成图，标记的区域为结合高分光学影像选取的实验样本，其中红色多边形表示倒塌建筑物，绿色多边形表示完好建筑物，蓝色则代表非建筑区。图5(b)为2010年5月6日玉树城区0.5 m分辨率光学遥感图像。图5(c)为基于光学图像解译的建筑物损毁程度参考结果。
3.1.1 非建筑区剔除
对极化数据Pauli分解并滤波处理，图6(a)是滤波后π/4偶次散射分量功率ω′。由图6(a)可见，非建筑区(道路、河流、裸地等)偏向于蓝色，而建筑区多呈现黄色或红色，建筑区的特征值整体上比非建筑区的值大。
根据第2.1节中描述的阈值确定方法确定区分建筑区和非建筑区的阈值ωt为–13.5(见图3中的黑色虚线)。剔除非建筑区的结果如图6(b)所示，其中黑色代表非建筑区，白色表示建筑区。本文方法不仅能剔除城区中东南向的宽河流，也较好地剔除了东西向较窄河流以及城区中心的道路，此外还能剔除部分城区中的裸地。
为了定量评价非建筑区剔除的效果，在建筑区和非建筑区分别选取999和987个像素来评估本文中非建筑区剔除的精度。建筑区和非建筑区的识别总精度为89.63%，其中非建筑区的检出率为83.18%，城区的检出率为96.00%。影响非建筑区剔除精度的主要原因是城区公园等区域由于π/4偶次散射分量较大而易被识别为建筑区。但最终待评估的建筑区不包含这种较大面积连续的区域，即这种情况对后续建筑物损毁评估不会造成大的影响。
3.1.2 倒塌建筑物提取
采用2.1节中描述的阈值确定方法确定损毁建筑物提取阈值，其中π/4偶次散射功率的variance纹理、contrast纹理，奇次散射功率的contrast纹理特征的阈值分别为20.34，40.00和55.70，如图4中的黑色虚线所示。分别用这3种纹理特征提取倒塌建筑物的结果如图7所示，三者均能较好地提取出主要的倒塌建筑，但也有一些存在差异的区域：绿虚线区表示提取结果相对准确的区域，红虚线区表示准确度相对偏低的区域。
为了分析综合3种纹理特征阈值划分对于完好建筑物和倒塌建筑物的区分能力，统计单纹理特征及综合三者后对样本点的区分精度，如表1所示。其中，综合纹理表示针对3个单独的纹理特征用投票法确定最终区分结果(即对区分结果的平均)。样本点在3个纹理特征组成的3维空间中的分布如图8所示。由表1和图8可知，识别的完好建筑物中不包含倒塌建筑物，识别的倒塌建筑物中包含极少量的完好建筑物。通过对3种纹理特征阈值划分的综合，对倒塌建筑物的检出率为100%，虚警率为7.06%；完好建筑物的检出率为92.40%，虚警率为0，能很好地区分完好建筑物和倒塌建筑物。综合3种纹理特征阈值划分后的区分能力明显优于单个纹理特征的阈值划分。
表 1 完好建筑物和倒塌建筑物的检出率
Table 1 The detection rate of intact and collapsed buildings
| /% | ||||
| 类型 | 检出率 | |||
| ω的variance | ω的contrast | u的contrast | 综合纹理 | |
| 倒塌建筑物 | 98.78 | 100.00 | 78.37 | 100.00 |
| 完好建筑物 | 85.65 | 80.17 | 93.67 | 92.40 |
3.1.3 建筑物损毁评估
由图7可见，3种纹理特征各有优势，为更准确地评估建筑物的损毁程度，综合利用这3种极化参数的纹理特征。
由式(6)基于图7(a)(b)(c)分别计算区域建筑损毁指数，三者取平均计算得到综合的区域建筑损毁指数。依据第2.3节中建筑物损毁程度的划分标准，将建筑物损毁程度划分为3个等级，得到本文方法的评估结果，如图9(d)所示。
对比实验H-α-ρ方法、总功率同质性纹理方法、综合NCCC和同质性纹理方法的结果如图9(a)-(c)所示。由于采用了相同的实验数据，图9(a)中H-α-ρ方法使用文献(郭华东 等，2010)中的阈值。根据第2.1节中描述的阈值确定方法，确定用总功率同质性纹理区分完好建筑物和倒塌建筑物的阈值为0.34，用NCCC识别倒塌建筑物的阈值为0.03。
对照参考图比较不同方法的评估结果可知，本文方法的建筑物损毁评估结果明显优于其他3种方法，更加接近参考结果。如图9(d)所示，本文方法不仅正确评估出了城区西部、西南部的严重损毁区，也基本准确地评估出了城区中部的轻度损毁建筑区；特别地，城区东北区方向大部分存在取向角的完好建筑的损毁程度得到正确评估。而H-α-ρ方法、总功率同质性纹理方法，将城区东北部的大片有取向角完好建筑区都错误地评估为重度损毁区(见图9(a)和9(b))；综合NCCC与同质性纹理的方法，能较好地评估该区域建筑物损毁程度(图9(c))。然而，对于城区南部存留有少数与方位向平行建筑物的严重损毁区，即图9(d)中绿虚线框区域，这3种方法存在不同程度的误评估；而本文方法能准确地评估出该严重损毁区。可见，本文方法受建筑物取向角的影响小，对于存在取向角的完好建筑物以及存留有少数与方位向平行的完好建筑物的损毁区，均能较准确地评估，表现出了更好的建筑物损毁评估能力。
进一步对建筑物损毁评估的结果进行定量评价，对不同损毁程度的建筑物从检出率、虚警率、总精度3个方面评价，不同方法评估结果的精度见表2。对于重度损毁建筑区，4种方法的检出率都较高，尤其是本文方法对重度损毁建筑的检出率为100%，同时虚警率也是4种方法中最低的。中、轻度损毁建筑区的评估精度相对较低，但本文方法的检出率仍然最高，虚警率也相对较低。从总精度而言，本文方法对于3种不同程度的建筑物损毁评估的效果最好，总精度达到74.39%。
表 2 不同方法的建筑物损毁评估结果精度评价
Table 2 Building damage assessment accuracy of different methods
| /% | ||||||||
| 方法 | 检出率 | 虚警率 | 总精度 | |||||
| 重度损毁 | 中度损毁 | 轻度损毁 | 重度损毁 | 中度损毁 | 轻度损毁 | |||
| H-a-p法 | 86.02 | 28.38 | 20.26 | 37.65 | 34.91 | 79.07 | 48.95 | |
| 同质性纹理法 | 73.26 | 17.77 | 24.50 | 49.39 | 76.51 | 62.31 | 42.55 | |
| 综合NCCC和同质性纹理法 | 75.32 | 15.32 | 51.45 | 37.81 | 59.14 | 57.56 | 51.98 | |
| 本文方法 | 100.00 | 44.55 | 66.15 | 15.62 | 40.37 | 31.67 | 74.39 | |
本文方法建筑物损毁评估的误差矩阵见表3，中、轻度损毁建筑区之间评估结果存在一定比例的混淆。由表3可见，本文方法评估结果中仍然存在少部分轻、中度损毁区被误评估为重度损毁区，如图9(d)中红虚线框区域。究其原因，一般被正确评估的完好建筑区主要为城区较繁华的地带的建筑物，房屋较为现代化，高低不一，如图10(a)所示，其方差和对比度纹理特征呈现出高特征值。而这些被误评的区域均为类似城边村的区域，房屋大都为老式房屋，低矮平整，排列紧凑，如图10(b)所示，其Pauli极化特征及其方差和对比度纹理特征与图10(c)所示的重度损毁建筑区较为相似，都呈现出低特征值，因此易被误评估。
表 3 本文方法建筑物损毁评估误差矩阵
Table 3 Building damage assessment matrix of the proposed method
| 评估\参考 | 重度损毁 | 中度损毁 | 轻度损毁 | 总计 |
| 重度损毁 | 10792 | 952 | 1046 | 12790 |
| 中度损毁 | 0 | 3006 | 2035 | 5041 |
| 轻度损毁 | 0 | 2790 | 6020 | 8810 |
| 总计 | 10792 | 6748 | 9101 | 26641 |
3.2 ALOS-1数据实验
为了进一步验证方法的适用性，利用2011年3月11日东日本大地震和海啸后的ALOS-1数据进行损毁建筑物评估。实验数据为2011年4月8日获取的日本宫城县石卷市的ALOS-1全极化数据，影像大小为494×245像素，入射角为23.83°，方位向进行8视处理。图11(a)为该实验区的PauliRGB合成图，图11(b)为根据高分光学数据解译的建筑物损毁程度图。
用本文方法对该区域建筑物损毁评估的结果如图11(c)所示。采用第2.1节中描述的阈值确定方法，确定π/4偶次散射功率的variance纹理、contrast纹理、奇次散射功率的contrast纹理特征区分倒塌建筑物和完好建筑物的阈值分别为10.00，16.50和15.00。选取图11(b)中标号为1—10的10个区块进行评估，其中1—4为重度损毁，5—7为中度损毁，8—10为轻度损毁。对轻度损毁、中度损毁、重度损毁的识别精度分别为100%、78.38%、65.01%，总精度为80.26%。
4 结 论
本文提出了一种采用震后单时相SAR数据、综合利用极化分解后多纹理特征的建筑物损毁评估方法。该方法利用Pauli分解的π/4偶次散射分量功率来区分非建筑区和建筑区，更完整剔除河流、道路和裸地等非建筑物区域，有效减少了非建筑区对建筑物损毁评估的影响。该方法综合利用Pauli分解后极化功率的三种灰度共生矩阵纹理特征：π/4偶次散射功率的方差、对比度纹理和奇次散射功率的对比度纹理，提出将极化特征与纹理特征嵌套融合以充分利用多类特征进行建筑物损毁程度评估。采用玉树RADARSAT-2数据和日本石卷ALOS-1数据的实验证明了本文方法对震后建筑物损毁评估的有效性，该方法对玉树城区和日本石卷城区损毁建筑评估的总体精度分别为74.39%和80.26%。与其他方法的对比实验表明，本文方法能降低取向角的影响，对大取向角的完好建筑区以及存留有少数与方位向平行的完好建筑物的倒塌区的评估更准确。当然，本文方法仍然存在一些不足：由于与倒塌建筑区纹理的相似性，一些低矮、紧凑的完好建筑物难以被正确评估；非建筑区剔除时，对城市公园等植被较多的非建筑区剔除效果不够理想。因此，如何提高对低矮、紧凑完好建筑物的评估效果、如何更准确剔除非建筑区以及进一步验证该方法对高空间分辨率极化SAR数据的适用性，是下一步的工作重点。
志 谢 感谢加拿大太空署与MDA公司的RADARSAT-2数据、日本JAXA的ALOS-1全极化数据以及德国宇航中心提供的玉树震后建筑损毁程度图。
参考文献(References)
-
Chen S W and Sato M. 2013. Tsunami damage investigation of built-Up areas using multitemporal spaceborne full polarimetric SAR images. IEEE Transactions on Geoscience and Remote Sensing, 51 (4): 1985–1997. [DOI: 10.1109/TGRS.2012.2210050]
-
Dekker R J. 2011. High-resolution radar damage assessment after the earthquake in Haiti on 12 January 2010. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 4 (4): 960–970. [DOI: 10.1109/JSTARS.2011.2164516]
-
Fan Y D, Wu W, Wang W, Liu M and Wen Q. 2016. Research progress of disaster remote sensing in China. Journal of Remote Sensing, 20 (5): 1170–1184. [DOI: 10.11834/jrs.20166171] ( 范一大, 吴玮, 王薇, 刘明, 温奇. 2016. 中国灾害遥感研究进展. 遥感学报, 20 (5): 1170–1184. [DOI: 10.11834/jrs.20166171] )
-
Gamba P, Dell’Acqua F and Lisini G. 2006. Change detection of multitemporal SAR data in urban areas combining feature-based and pixel-based techniques. IEEE Transactions on Geoscience and Remote Sensing, 44 (10): 2820–2827. [DOI: 10.1109/TGRS.2006.879498]
-
Guo H D, Li X W and Zhang L. 2009. Study of detecting method with advanced airborne and spaceborne synthetic aperture radar data for collapsed urban buildings from the Wenchuan earthquake. Journal of Applied Remote Sensing, 3 (1): 031695 [DOI: 10.1117/1.3153902]
-
Guo H D, Wang X Y, Li X W, Liu G, Zhang L and Yan S Y. 2010. Yushu earthquake synergic analysis study using multi-model SAR datasets. Chinese Science Bulletin, 55 (13): 1195–1199. [DOI: 10.1007/s11434-010-0292-2] ( 郭华东, 王心源, 李新武, 刘广, 张露, 阎世勇. 2010. 多模式SAR玉树地震协同分析. 科学通报, 55 (13): 1195–1199. [DOI: 10.1007/s11434-010-0292-2] )
-
Hao H M, Zhang Y H, Shi H Y and Huang J B. 2012. Application of test statistic method in fully polarimtric SAR change detection. Journal of Remote Sensing, 16 (3): 520–532. [DOI: 10.11834/jrs.20120384] ( 郝洪美, 张永红, 石海燕, 黄金波. 2012. 统计假设检验方法在全极化SAR变化检测中的应用. 遥感学报, 16 (3): 520–532. [DOI: 10.11834/jrs.20120384] )
-
Haralick R M, Shanmugam K and Dinstein I. 1973. Textural features for image classification. IEEE Transactions on Systems, Man, and Cybernetics, SMC-3 (6): 610–621. [DOI: 10.1109/TSMC.1973.4309314]
-
Lee J S and Pottier E. 2009. Polarimetric Radar Imaging: From Basics to Applications. Boca Raton: CRC Press Taylor and Francis Group: 214–215
-
Matsuoka M and Yamazaki F. 2004. Use of satellite SAR intensity imagery for detecting building areas damaged due to earthquakes. Earthquake Spectra, 20 (3): 975–994. [DOI: 10.1193/1.1774182]
-
Shi L, Sun W D, Yang J, Li P X and Lu L J. 2015. Building collapse assessment by the use of postearthquake Chinese VHR airborne SAR. IEEE Geoscience and Remote Sensing Letters, 12 (10): 2021–2025. [DOI: 10.1109/LGRS.2015.2443018]
-
Yonezawa C and Takeuchi S. 2001. Decorrelation of SAR data by urban damages caused by the 1995 Hyogoken-Nanbu earthquake. International Journal of Remote Sensing, 22 (8): 1585–1600. [DOI: 10.1080/01431160118187]
-
Zhai W, Shen H F, Huang C L, Huang C L and Pei W S. 2016. Building earthquake damage information extraction from a single post-earthquake PolSAR image. Remote Sensing, 8 (3): 171 [DOI: 10.3390/rs8030171]
-
Zhai W and Huang C L. 2016. Fast building damage mapping using a single post-earthquake PolSAR image: a case study of the 2010 Yushu earthquake. Earth, Planets and Space, 68 : 86 [DOI: 10.1186/s40623-016-0469-2]
-
Zhang H Z, Wang Q, Zeng Q M and Jiao J. 2015. A novel approach to building collapse detection from post-seismic polarimetric SAR imagery by using optimization of polarimetric contrast enhancement//Proceedings of 2015 IEEE International Geoscience and Remote Sensing Symposium. Milan, Italy: IEEE: 3270–3273 [DOI: 10.1109/IGARSS.2015.7326516]
-
Zhang J F, Xie L L and Tao X X. 2002. Change detection of remote sensing image for earthquake damaged buildings and its application in seismic disaster assessment. Journal of Natural Disasters, 11 (2): 59–64. [DOI: 10.3969/j.issn.1004-4574.2002.02.010] ( 张景发, 谢礼立, 陶夏新. 2002. 建筑物震害遥感图像的变化检测与震害评估. 自然灾害学报, 11 (2): 59–64. [DOI: 10.3969/j.issn.1004-4574.2002.02.010] )
-
Zhao L L, Yang J, Li P X, Zhang L P, Shi L and Lang F K. 2013. Damage assessment in urban areas using post-earthquake airborne PolSAR imagery. International Journal of Remote Sensing, 34 (24): 8952–8966. [DOI: 10.1080/01431161.2013.860566]
