基于激光测高卫星全波形数据的地表覆盖分类应用——以北京市城区为例

权学烽，唐新明，李国元，高显连

(1.兰州交通大学，兰州 730070；2.自然资源部国土卫星遥感应用中心，北京 100048；3.甘肃省地理国情监测工程实验室，兰州 730070；4.国家林业局卫星林业应用中心，北京 100714)

0 引言

1971年美国的阿波罗-15号搭载了第一台星载激光测高仪，并在随后的时间里大力发展卫星激光测高技术[1-7]。2003年美国发射了首颗对地观测激光测高卫星，冰、云和陆地高程卫星ICESat/GLAS(ICE，Cloud and land Elevation Satellite/Geoscience Laser Altimeter System)，该卫星于2009年停止数据接收，在这期间共接收了近2亿个激光测高点，为土地覆盖分类研究提供了丰富的数据。近几年来，随着对地观测技术的发展和国家战略的需求，国家高分辨率对地观测系统重大专项以及国家民用空间基础设施中均明确提出我国将发射高分七号(GF-7)卫星和陆地生态系统碳监测卫星[8]。这些卫星都将搭载激光测高仪，为我国的星载激光测高数据应用研究提供帮助。2016年5月30日，资源三号02星顺利发射，02星上搭载了国内首台对地观测的激光测高仪，主要用于测试激光测高仪在对地测量时的精度，以及探索激光测高数据与光学遥感影像结合实现无控立体测图的可行性[9]。目前，国外围绕基于星载全波形激光测高数据在土地覆盖分类应用开展了较多的研究，并且以美国的ICESat数据为主，国内仅有少数团队在进行星载激光测高数据的应用研究，其中森林生物量反演和湖泊监测等方向的研究人员较多，而土地覆盖分类方向的研究还处于刚刚起步的阶段。

土地覆盖分类应用研究主要目的是通过激光测高仪接收到地表回波的波形信息来进行地表覆盖物的分类，提升分类准确率一直是该研究领域的热点与难点。目前该方向的研究主要通过两个部分来提升分类准确率，一部分是在波形处理时使返回波波形更为精确以及提高回波波形的利用率，挖掘有效的波形参数，其中波形处理的方法有利用高斯分解提取波形特征参数的波形特征参数法和将波形转化为累计分布函数(cumulative distribution function，CDF)曲线的曲线匹配法2种方法；另一部分是通过对多种分类器的研究与改进来提高分类的准确率，目前应用于星载激光测高数据在地表覆盖分类的主要分类方法有K均值法(K-mean，KM)、围绕中心剖分法(partitioning around medoids，PAM)、模糊C均值法(fuzzy C-means，FCM)、决策树算法(decision tree，DT)、最大似然分类(maximum likelihood，ML)、人工神经网络(artificial neural networks，ANN)、随机森林(random forest，RF)以及支持向量机(support vector machines，SVM)。Duong等人通过对原始波形进行波形预处理后利用波形特征参数法进行了土地覆盖分类，实验分类准确率为73%，表现出良好的一致性[10]。大多数国外学者在进行土地覆盖分类研究时也使用该方法，通过对不同地区的实验发现，选取适当的分类器可以提高分类的准确率，但是当分类的种类增多时，分类的准确率出现下降，且该方法是对连续的波形进行参数提取变为离散的数据，在面对多种相似但不同的种类时，该方法可能损失了部分有用的信息导致分类的准确率下降，但是该方法在大多数情况下都是适用的[11-15]。Zhou等人提出了基于KS距离分类的CDF曲线匹配法，实验结果表明基于KS距离分类的土地覆盖分类准确率可以达到87.2%，Kappa系数为0.80[16]。这种方法优于使用高斯分解参数的土地覆盖分类方法3.5%，该方法在一定的区域内提供了比特征参数法更高的分类准确率，这可能是由于该方法是对完整波形的使用，获得了更多的分类信息，从而使分类准确率上升。2016年Zhou等人利用全波形数据来测试曲线匹配法区分不同垂直结构物体的能力，在基于原来的曲线匹配算法的基础上扩展了2组新的曲线匹配方法，分类准确率有所上升[17]。

本文使用ICESat/GLAS在北京地区的全波形数据，采用改进后的CDF和KS距离算法，开展土地覆盖分类实验研究，以验证星载全波形激光数据在土地覆盖分类的应用前景，为我国后续的激光测高卫星在土地覆盖分类应用提供参考。

1 实验区和实验数据

1.1 实验区概况

实验区位于中国北京市市区周边区域(115.4°E～117.5°E，39.4°N～41.1°N)。北京市中心位于116°23′E、38°54′N，位于华北平原西北边缘，与天津相邻且被河北省环绕包围，如图1所示。

图1 北京市GLAS足印点分布图

由于全波形数据对于复杂地物分类存在一定局限性，相似地物类型在分类时的效果较差[11]，因此本研究将实验区地表类型分为4类，分别为裸地、水体、树木和建筑物。该区域内的建筑物构造复杂，既有较高的楼房也有低矮的棚户型建筑，因此分类难度较高。由于北京周边地形复杂、有较高的高程变化，为了减弱地形对于分类结果的影响，故选取北京市城区进行土地覆盖分类实验。

1.2 实验数据

1)ICESat/GLAS数据。实验中使用ICESat/GLAS数据为分类数据，ICESat/GLAS是首颗对地观测激光测高卫星，且大量国外学者已经验证了该数据的可靠性，因此使用该数据进行土地覆盖分类得到的结果具有一定的可靠性。本实验使用北京地区的2003—2008年的GLA01、GLA05和GLA14。GLA01数据产品包含发射和接收的波形数据、光斑号、时间、经纬度、脉冲能量等信息；GLA05数据产品包含波形特征参数、激光足印的长短轴参数和其他计算表面倾斜等地形特征时必须用到的参数；GLA14数据产品包含经波形分解得到的高斯波信息以及地表的高度信息、经纬度和高程信息。

2)验证数据GLC30。本实验使用2003—2008年Landsat-7历史遥感影像及GlobeLand30数据进行验证，通过对足印点与遥感影像的叠加分析可以精确地判读该区域内的土地覆盖信息。GlobeLand30分类数据利用的影像为30 m多光谱影像，包括美国陆地资源卫星(Landsat)TM5、ETM+多光谱影像和中国环境减灾卫星(HJ-1)多光谱影像。除了多光谱影像外，还使用了大量的辅助数据和参考资料，以支持样本选取、辅助分类等工作。在GLC30数据集中代码10为耕地、20为森林、60为水系、80为人造覆盖物、90为裸地。本实验将耕地与裸地合并为裸地类，因为其波形数据较为相似，且当耕地农闲时其地表特性与裸地相似；森林和部分单独的树木分类为树木类。

2 分类方法

本实验为3组实验，第1组实验为初步实验，验证CDF曲线KS距离匹配法在土地覆盖分类方向应用的可行性，并显示不剔除坡度、地形因素的分类准确率；第2组实验为剔除坡度与地形因素后的分类实验，与第1组实验组成对比实验，显示出坡度与地形对地表覆盖分类准确率的影响；第3组为本文改进方法的实验，该实验在CDF曲线与KS距离匹配法的基础上加入了波峰数和返回能量2个特征值，期望能进一步提升分类的准确率。本文改进方法实验流程如图2所示。

图2 本研究方法流程图

2.1 数据预处理

GLAS数据预处理是数据应用的前提和保障，因此需要以下几个步骤对数据进行处理校正后才能进行进一步的应用。

①电压值转换。GLAS数据被压缩为1～256的计数波形，电压值转换是将计数波形转换为真实波形的过程。

②波形正则化。在对比分析不同脉冲和环境下的波形时需要利用波形正则化将其统一。

③数据解压缩。GLAS数据是将连续的1 000帧信号压缩为544个(陆地)或200个(海洋)采样值，数据解压缩过程就是该过程的逆过程。

④波形滤波。波形滤波就是将波形中存在的噪声进行去除。

在进行完上述标准步骤后，对足印点进行进一步的筛选，首先将返回能量值小于0.1的数据进行剔除，这是因为该数据受到云层的影响较大，之后对坡度大于10°的数据进行剔除，以消除坡度对地物分类的影响[18]。

2.2 CDF转换

CDF转换是指将经过数据预处理的波形数据转化为累计分布函数图像的过程，可以由公式(1)进行转换。

(1)

式中：jmin为阈值设置后波形开始的值(即开始接受能量的时刻)。

2.3 KS距离匹配

KS距离匹配(Kolmogorov-Smirnov distance text，KS)是用来测试样本曲线与实验曲线之间的相似程度，通过实验曲线与样本曲线之间的对比得到实验曲线与更接近于哪一组样本曲线。KS距离匹配方法由公式(2)、公式(3)表示。

(2)

(3)

2.4 返回能量对比

通过对不同足印点反射能量大小进行对比可以进一步准确地对地物进行分类，有效地区分裸地与水体、建筑物与树木。返回能量由公式(5)近似表达。

(4)

式中：vi为ti时刻的电压值，其中t为计数间隔，不影响不同地物返回能量之间的对比，所以可化简为

(5)

通过对已经分类后的各个种类进行能量计算后，对每个种类计算其期望值作为能量阈值，按大小排列。当已知A类的能量阈值a，B类的能量阈值b，当aB时，将a1的数据分为B类。当存在3类时，C类的能量阈值为c且ab且a1c时将a1分为C类。

3 实验与分析

3.1 样本训练实验

使用目视解译法在筛选后的数据中选取20个具有良好代表性的激光点数据(每种类型5个点)之后，随机提取200个激光点作为测试数据进行分类，使用KS距离匹配法进行分类，分类结果如表1所示。

表1 样本分类结果

其中有29个建筑物为正确分类，84个裸地为正确分类，4个水体为正确分类，22个植被为正确分类，将这139个波形数据作为样本数据进行之后的实验。

图3为具有代表性的典型地物的地面真实影像与波形图象、CDF曲线的匹配图。

(图中红色椭圆为近似激光足印点)图3 影像、波形、CDF匹配图

3.3 CDF方法

本实验使用曲线匹配法进行分类，实验数据选取了2003—2009年的GLAS点共8 600个点进行实验。该数据集离散地散落于整个北京市，整个数据集之间的空间关联性较差，其分类结果如表2所示。

表2 实验1分类结果

其准确率为0.67，其Kappa系数为0.415。通过对表2中错误分类点的分析发现，当选取的实验点数据未处于平坦地形且空间关联性较差时，较高地形的树木会错误地分类为建筑物，相对应的较低地形的建筑物也会错误地分类为树木，辨别建筑物与树木主要是因为建筑物的CDF曲线上升较快与树木CDF曲线不同，且建筑物与树木的起波点不同且转折点位置不同，但是由于地形原因导致树木与建筑物的起波点相似，使得其整体不同差距变小，而不同地形的建筑物在起波点的差距导致其具有较大的不同性，使得树木错误地分类为建筑物，建筑物错误地分类为树木。实验发现裸地和水体之间，裸地和建筑物、树木之间也有错误分类的情况，分析后发现这可能是由于建筑物或树木的地表反射能量较低，形成的上升曲线与裸地或水体的曲线较为相似从而错误地分类，且也可能是由于裸地的起波点与树木和建筑的起波点相似造成的错误分类。而当坡度较大时会造成无法判断激光点的最后一个波峰是否为地面反射的波峰的情况，而且也会造成波形展宽的现象。

通过以上分析可以发现在不同地物之间的分类其起波点位置在分类中显得尤为关键，也就是说地形和地物高程的不同较为关键，因此在进行地物分类时要选取同一地势、坡度较小的区域进行分类。

3.4 去除高程影响的CDF方法

该实验继续使用与实验1相同的方法，但在激光点选取时选取了同一地形、坡度不大于10°的激光点数据[18]，并且对范围内的激光点数据进行了进一步筛选，最后得到筛选后的2 000个激光足印点，分类结果见表3所示。

表3 实验2分类结果

其准确率为0.81，其Kappa系数为0.58，分类准确率和Kappa系数有较高的提升。因此，坡度和地形对于分类准确率还是有较大的影响。本实验还通过对被筛选掉的数据分析发现，大多数处于高坡度的地表覆盖类型都是裸地或者树木，基本没有建筑物，即使存在建筑物也为临时搭盖的棚户房。

通过对错误分类的数据进行分析发现，在同一地形下，建筑物和树木出现相互错误分类的情况，通过在高分辨率遥感影像上发现，这是由于部分建筑物为低矮建筑物，只有2～4 m，与树木高度一致，导致起波点相似，在分类时不能很好地判别。水体、裸地与建筑物、植被之间也出现少量的分类错误现象，这可能也是由于建筑物或树木的高度较低，且地表返回能量较低导致形成的曲线与裸地、水体的CDF曲线较为相似导致的。

3.5 本文方法

实验在Zhou的方法的基础上加入了波峰数与返回能量值的参数分类。首先对CDF曲线进行分类，分为波峰(水体、裸地)与多波峰(建筑物、树木)两类，该过程也可以在CDF曲线转换之前进行。在进行KS距离匹配法之后，进行返回能量判别。其分类结果如表4所示。

表4 本文方法分类结果

准确率为0.90，Kappa系数为0.767，表现出较好的分类准确率，但是Kappa系数较低。该结果验证了该改进方法的可行性。可以发现多数错误的分类已被分为正确的类型。在本实验中发现当分类数较少时，其返回能量阈值之间的间距较大，不会产生重叠的情况，当分类种类较多时，返回能量阈值之间的间距会变小，区别性会降低，该方法的分类效果将会变低，但是分类数与分类准确率之间的关系还需进一步实验研究。实验结果中仍有个别的类型还是处于错误分类之中，通过对比同时相的影像数据发现，其中建筑物错误地分类为裸地是由于验证样本的时间与实验数据之间的不一致性导致的。通过对比历史影像和GLC30分类图发现，GLC30为2010年的数据生成的地物类别地图，而北京市在2003—2010年之间发生了较大的变化，其中北京南部郊区等地的地表覆盖类型都由裸地转变为建筑物，而使用GLAS数据可能正确地分类裸地但是在判断时错误地判断为建筑物。除去这部分因素之后分类准确率还会进一步上升，其他的错误分类还需进一步研究。

通过表5可以清晰地发现去除坡度和地形因素地表覆盖分类准确率的上升，改进后的方法在准确率和kappa系数上都有较大的提升，表明了该方法的有效性。

表5 实验结果对比

4 结束语

本研究采用GLAS数据，将CDF曲线匹配法与KS距离匹配法相结合的方法，验证了卫星激光测高全波形数据在地表覆盖类型分类中的可应用性，并对该方法进行了适当的改进以提升其分类准确率。实验结果表明该方法可以有效提高分类的准确率，但是还存在一定的问题。首先该方法去除了坡度与地形的影响，这说明该方法仅适用于平原甚至是城市地区，在坡度较大的地区无法进行有效的分类，对分类的范围有了制约；第二，在进行阈值设定时，需要分类类型之间具有较大的区别，因为相似类型的覆盖物的能量阈值极为相似，不能有效进行二次判别；第三，本实验选取的建筑物与植被阈值的可变性。通过实验发现，就树木与建筑物的返回能量而言，树木的返回能量一般大于建筑物的返回能量，但是这与树木冠层的大小、树木的高度有关，而区域性的树木和城区范围内的树木一般为同批次同种类栽种，其冠层大小与高度相似，因此其能量阈值设定较为规律与简单。而当分类区域为野外时，树木生长的时间不一致并且种类混杂，其能量阈值的选取较为困难，有待进一步实验进行论证。