徐宏根，宋 妍

(1.武汉地质调查中心，武汉 430205;2.中国地质大学(武汉)信息工程学院，武汉 430072)

0 引言

遥感图像变化检测是近10 a以来遥感应用研究的热点问题之一，在军事目标侦查、灾害评估、GIS数据库更新等领域有着广阔的应用前景。本文以覆盖城市的高空间分辨率遥感图像为研究对象，旨在解决城市场景内由阴影所引起的误检测问题。

已有许多学者研究过如何降低遥感图像中阴影对变化检测的影响，例如:Shintaro[1]在同名核线上搜索房顶的信息，并根据阴影模型去除因阴影造成的伪变化;季顺平[2]利用高斯背景模型提取阴影，并在对图像完成阴影补偿后再进行变化检测。对于上述方法，阴影提取的准确性是关键问题。对这个问题，王树根[3]认为阴影的提取可分为低级处理与中级处理2步:首先依据光谱特征提取阴影的初始区域;然后依据阴影的几何特点对其进行确认，并提出用“基于RGB模型空间的阴影检测”方法检测彩色航空图像中的阴影。针对高分辨率遥感图像中阴影面积较大的特点，郭海涛等[4]还提出从光谱和几何2个方面进行约束和提取阴影的方法。

针对2期遥感图像之间由于阴影的变化所造成的伪变化，本文采取阴影提取后对初始检测结果进行阴影补偿的方法来提高变化检测的精度。该方法的核心是阴影提取。在提取2期遥感图像的阴影后，即可根据2期图像的阴影差分结果得到2期图像因阴影造成的变化，从而去除初始检测结果中阴影所造成的影响。在阴影提取时，本文借鉴面向对象的分类思想，采用先分割、后提取的方法完成阴影检测。最后通过实验证明本文方法的有效性和精度。

1 阴影检测

阴影检测是剔除阴影对变化检测结果负面影响的先决条件。下面首先分析阴影形成原因，并总结阴影的光谱和几何方面的性质。

图像中的阴影是成像光源被目标物体完全或部分遮挡所形成的，只要传感器与光源所在位置(或方向)不重合，就会产生阴影[4]。从光谱特征而言，阴影的光谱值较低，根据这一特点可将阴影与其他色调较亮的物体区分开来;其次，由于天空中大气介质对蓝光分量有较强的散射作用，蓝波段的数据在阴影区域的方差较大，依此也可区分阴影与其他一些光谱值较低的物体。另外，植被与阴影在可见光波段均呈暗色调，难以通过可见光谱段的光谱进行区分;但在近红外波段植被呈现亮色调，而阴影依旧保持暗色调，因此利用近红外波段的数据可以较好地区分植被与阴影。

根据阴影的上述性质，本文采用以下步骤提取阴影:①分别对2期遥感图像进行波段选择，选择出适合阴影检测的波段;②运用图像分割的方法，对步骤①所选择的波段数据进行分割;③根据各个分割块的均值将各个分割块分为3类:阴影、疑似阴影和非阴影;④进一步根据分割块在蓝波段的方差确定最终的阴影。

1.1 波段选择

为提取阴影，需要从遥感图像的原始波段中选择出阴影可分性较好的波段。鉴于近红外波段(nearinfrared，N)与蓝波段(blue，B)的数据对提取阴影有较好的作用，在选择阴影可分性较好的波段的基础上，本文继续在红波段(red，R)与绿波段(green，G)中进行选择。运用遥感专业软件Erdas分别按照(B+R+N)、(B+G+N)以及(B+G+R+N)等3 种方式进行波段合成，并计算3幅合成图像中阴影与其他地物类型的转换散度。Swain等[5]指出:“若2类别间的转换散度大于1 900，则这2类可以被很好地分开;若转换散度在1 700～1 900之间，则这2类是可分的，分割效果尚可;若转换散度小于1 700，则这2类的可分性很差。”因此，从波段数目与转换散度可分性2个方面综合考虑，选择出适宜的假彩色合成图像进行后续分析。

1.2 图像分割

对于高分辨率遥感图像而言，传统的面向像元的分类方法会使分类结果中产生“椒盐效应”，导致分类精度低，而面向对象的图像分类方法是解决上述问题的有效方法。但是，面向对象的分类方法对图像分割的结果有较强的依赖。在彩色航空图像的阴影检测中，已有人探索将均值平移(mean shift)算法用于图像分割，取得了比传统方法更好的效果[6]。因此，本文也采用均值平移算法来实现对原始图像的分割。

均值平移算法实际属于非参数估计密度函数的方法，最早由 Fukunaga[7]在 1975 年提出;Comaniciu[8]则证明了均值平移算法具有严格的收敛性，并且将均值平移方法运用到图像分割中。均值平移算法将图像分割转化为寻求图像特征空间模式的过程，并分为概率密度估计、模式搜索、图像滤波和图像分割4个连续的步骤(其具体的计算流程见文献[8])。

1.3 阴影提取

完成图像分割后，即可依据分割块的光谱和几何性质提取阴影区域。笔者认为，原始遥感图像的光谱特征符合高斯混合模型(Gaussian mixture model，GMM)，因此可依据贝叶斯最小错误概率准则、根据各分割块的均值提取初始的阴影区域。但由于阴影与图像中的暗色调地物容易混淆，故需要进一步计算各分割块在蓝波段的方差，从而对初始的阴影提取结果进一步约束，以得到真正的阴影区域。

1.3.1 运用分割块均值提取初始阴影

为了利用各个分割块的光谱信息提取初始的阴影区域，首先采用GMM对整景遥感图像建模。整景图像可表示为

式中:zi为第i个高斯分布成分的权重因子，且满足归一化条件

和zi＞0;r为图像所包含的地物类型的数目;fi(x)为第i个高斯成分的密度函数，即

其中:mi为第i个成分的均值向量;σi为第i个成分的协方差矩阵;d为该成分的特征空间维数。

采用遗传K－均值算法、结合期望最大化(expectation－maximization，EM)算法求解出各个类别的统计参数[9－10]，即

然后统计各分割块的均值向量mn(1≤n≤N，N为分割块的总数)，依据Bayes公式计算分割块的均值向量属于不同类别的后验概率，并根据最小错误概率准则判断该分割块的类属，从而完成阴影的初始提取。

1.3.2 运用分割块方差约束初始阴影

在上述阴影提取的初始结果中，不仅包含了阴影，而且包含了水体等暗色调地物。为了去除干扰、得到真正的阴影，计算步骤①中各分割块在蓝波段中的方差，以蓝波段方差值作为各分割块的特征。阴影区域的方差较大，呈现较亮的色调;而水体则呈现较暗的色调。将各分割块的方差值依照从小到大的顺序排列，并人工设置阈值t，则所有方差大于t的分割块为阴影，而方差小于t的分割块为“非阴影”。阈值t的设置要通过反复试验的方法确定，以便尽可能地剔除非阴影和水体。

1.4 降低误检测

在提取出阴影区域后，将2期遥感图像的阴影提取结果直接相减，得到“阴影差分图”。该图反映了由于2期图像中阴影的变化所带来的伪变化，是导致城区遥感图像变化检测精度低的主要原因之一。为了去除这种伪变化，将阴影差分的结果与原始的变化检测结果作逻辑运算:如果一个像元在初始变化检测时被判断为“变化”的像元，同时在阴影差分图上也属于变化像元，则该像元并不是真变化像元，应将该像元改变为未变化像元;只有当像元在阴影差分图上属于未变化像元时，才保留初始变化检测的判断。

2 实验与结果

选择2003年和2004年获取的深圳地区IKONOS多光谱图像对(图1(a)和(b))作为实验数据，2景图像在进行变化检测之前经过相对辐射校正和几何配准，配准误差小于0.5像元。图像中的主要变化类型为城区内的植被变为建设用地，并有建设用地变为建筑。选择场景内变化的像元(图1(c)中的1 466个红色像元)和未变化的像元(图1(c)中的2 583个绿色像元)用于计算变化检测精度。图1(a)和(b)中蓝色椭圆内为城区内典型高楼及阴影区，在后续实验中将用于对比阴影去除前、后的变化检测效果。

图1 深圳地区IKONOS多光谱图像与差分图像Fig.1 IKONOS multi－spectral images and difference image in Shenzhen district

图2 基于灰度信息的图像变化检测结果Fig.2 Change detection result from gray value in the image

图2 是根据2期图像的光谱差分信息，采用顾及上下文信息的模糊聚类方法进行变化检测的结果[11]。该方法具有无需估计参数、自动化程度高、效率高等优点，变化检测总体精度为87.68%，Kappa系数为 0.72。

由实验结果可以看出，图像变化检测结果中存在着由于阴影的变化而造成的误检测，因而降低了变化检测的精度。具体地看，图3是对图1的IKONOS图像中蓝色椭圆区域裁剪后放大3倍的结果(图3(a)为前期图像，(b)为后期图像，(c)为仅依据灰度信息的变化检测结果)，用白色椭圆标出的2个区域被判断为变化区域。然而，这种变化是由阴影的变化造成的，应该予以剔除。

图3 由阴影引起的误检测Fig.3 Error detection caused by shadow

2.1 波段选择结果

散度是一种基于类别概率密度函数的可分性判据，设有类别 i，j，则其类间散度为

式中:si和mi分别为第i类的协方差矩阵和均值矢量;sj和mj分别为第j类的协方差矩阵和均值矢量。实际运算中，常采用变换散度TDij代替式(5)，即

对实验图像对进行判读，2期图像中所包含的典型地物有阴影、建筑物、道路、水体和植被。在3种波段合成的结果中，分别计算出阴影与其余4类地物的变换散度(表1)。

表1 不同波段组合的地物可分性比较Tab.1 Comparison among different surface feature separabilities of different band combinations

Swain[5]指出，“若 2 类的变换散度大于 1 900，则这2类可以被很好地分开;若在1 700～1 900之间，则这2类的可分性尚可;若小于1 700，则这2类的可分性将很差。”

由表1可以看出，对前期图像(图3(a))而言，3种波段合成方法都可以较好地将阴影与其余地物区分开，其中(B+R+N)波段组合在波段数目较少的情况下与(B+G+R+N)4个波段的组合有着相似的可分性结果。对后期图像(图3(b))而言，因为图像中云量较大，使得阴影与色调较暗的水体较难区分;从波段合成后的地物转换散度可分性结果看，(B+G+N)波段组合的波段数目较少，与4个波段合成结果的地物可分性相差无几。因此，下面的实验对前期图像采用B，R，N这3个波段的数据进行合成，而后期图像则采用B，G，N这3个波段的数据进行合成。

2.2 阴影提取结果

运用均值平移算法对2期图像进行分割，并依据GMM模型对分割后图像建模，根据均值信息提取初始的阴影(图4(a)和(b));依据阴影提取的初始结果，计算各个初始阴影块中蓝波段的方差，得到方差图像(图4(c)和(d));手工设置阈值t，认为所有方差大于t的分割块为阴影(图4(e)和(f))，而方差小于t的分割块为“非阴影”。

图4 基于光谱特征的阴影初始提取结果Fig.4 Initial extraction results based on spectral feature

2.3 阴影补偿变化检测结果

运用阴影差分方法补偿遥感图像初始变化检测的结果见图5。

图5 阴影补偿变化检测结果Fig.5 Change detection results compensated by shadow

图5 (a)是2期图像中阴影的差分结果，其中白色像元代表因阴影而产生的变化;图5(b)是运用阴影补偿图像初始变化检测的结果;为了更清晰地显示补偿效果，将图1(a)和(b)中蓝色椭圆的区域放大3倍，即为图5(c)和(d);图5(e)为补偿后的变化检测最终结果，其中白色像元为运用顾及上下文信息的模糊C均值方法检测时初始判断为变化的像元，黑色像元为该方法判断为未变化的像元，灰色像元为原始判断为变化、而由阴影差分判断为“伪变化”的像元。从图5中可以明显地看到由阴影造成的变化已被补偿掉。根据实验区的数据统计，阴影补偿后结果的总体精度为91.41%，Kappa系数为0.81;阴影补偿后总体精度比前述模糊聚类方法提高了3.73%，Kappa系数提高了0.09。

3 结论

利用高分辨率遥感图像在城区场景内进行变化检测，要重点解决由于阴影的存在引起的误检测。本文提出的运用阴影的光谱信息和几何信息、结合面向对象分类方法提取阴影和对初始的变化检测结果进行阴影补偿的方法，可以较好地去除因阴影引起的变化误检测，有效地提高高分辨率遥感图像变化检测精度。

[1]Shintaro W.A method for change detection of buildings using epipolar constraint from aerial images taken at different positions[J].Joho Shori Gakkai Kenkyu Hokoku，2000(50):33－40.

[2]季顺平.基于多时相航空影像的建筑物变化检测方法研究[D].武汉:武汉大学，2007.Ji S P.Building change detection research based on multi－ temporal aerial images[D].Wuhan:Wuhan University，2007.

[3]王树根.正射影像上阴影和遮蔽的成像机理与信息处理方法[D].武汉:武汉大学，2003.Wang S G.Imaging mechanism and information processing approach of shadows and occlusion on orthophotos[D].Wuhan:Wuhan University，2003.

[4]郭海涛，徐 青，张保明.多重约束下的建筑物阴影提取[J].武汉大学学报:信息科学版，2005，30(12):1059－1062.Guo H T，Xu Q，Zhang B M.Shadow extraction of building based on multiple constraints[J].Geomatics and Information Science of Wuhan University，2005，30(12):1059－1062.

[5]Swain P H，Shirley M D.Remote sensing:The quantitative approach[M].New York:McGraw Hill Book Company，1978.

[6]王军利.基于彩色航空影像的阴影检测算法研究[D].武汉:武汉大学，2005.Wang J L.The research of shadow detection algorithms based on color aerial images[D].WuHan:Wuhan university，2005.

[7]Fukunaga K.The estimation of the gradient of a density function with application in pattern recognition[J].IEEE Transaction of Information Theory，1975，21(1):32－40.

[8]Comaniciu D.Mean shift:A robust approach toward feature space analysis[J].IEEE Transactions on Pattern Analysis and Machine Intelligence，2002，24(5):603－619.

[9]Dempster A P.Maximum likelihood from incomplete data via the EM algorithm[J].Journal of the Royal Statistical Society，1977，39(1)－38.

[10]Krishna K，Narasimha M M.Genetic K－ means algorithm[J].IEEE Transactions on Systems，Man，and Cybernetics，Part B，1999，29(3):433－439.

[11]张 路，廖明生.一种顾及上下文的遥感影像模糊聚类[J].遥感学报，2006，10(1):58－65.Zhang L，Liao M S.Contextual fuzzy clustering of remote sensing imagery[J].Journal of Remote Sensing，2006，10(1):58－65.