面向遥感影像智能分类的海量样本数据采集方法

程 滔，吴 芸，郑新燕，杨 刚，白 驹

(1.国家基础地理信息中心，北京 100830；2.中国地图出版社，北京 100045)

随着计算机硬件的发展与计算能力的提升，大数据计算与分析、人工智能、深度学习等先进技术取得了实质性的快速进步[1-2]，并在商品识别、人脸识别、目标检测、敏感图识别、人类活动分析等实际应用场景中发挥了重要作用[3]，应用潜力十分巨大。

基于深度学习等技术的遥感影像智能分类是近年来遥感领域研究的热点[4-6]。这些技术的实质，是通过构建具有很多隐层的机器学习模型和海量的训练数据来学习更有用的特征，从而提升分类准确性。因此，遥感影像智能分类需要大规模训练样本和测试样本数据，并且样本数据准确度越高，通过训练进行识别和分类的精度就越高[7]。这里的样本数据指遥感影像和遥感影像对应的明确的地物目标组成的一一对应的关系数据组。

地理国情监测重大专项形成的覆盖全国范围的高分辨率遥感影像、高精度地表覆盖分类产品[8]，为遥感影像智能分类海量样本数据采集提供了丰富、可行的数据资源。基于遥感影像数据与地表覆盖分类栅格化数据，利用空间位置匹配技术，构建遥感影像与地表覆盖分类标签数据对，通过大规模采集，可为遥感影像智能分类提供精确的样本数据[9-10]。

本文针对遥感影像智能分类对高质量、大规模样本数据的需求，开展海量样本数据采集方法研究，解决采集过程中的各项技术难点；同时开发样本数据自动采集软件。利用整套方法，完成全国尺度海量样本数据的采集，并基于多个不同区域的采集过程与样本数据成果，对方法的实用性、运算性能进行评估。

1 数据源

1.1 高分辨率遥感影像

地理国情监测使用的遥感影像数据空间分辨率标准为全国优于2.5 m、重点区域优于1 m，以资源三号测绘卫星获取的影像为主，其他国产高分辨率遥感影像(包括高分一号、高分二号、北京二号、天绘一号等)及商业卫星遥感影像(包括WorldView-1/2、GeoEye、QuickBird、pléiade-1A/1B、SPOT 6/7等)进行补充，部分地区根据需要开展航空摄影、无人机航空摄影获取影像。原始影像各波段位深为16 bit，分幅数字正射影像各波段位深为8 bit[11]。

地理国情监测的标准时点为6月30日，这种多源遥感影像协同满足地理国情监测需求的现状，使得遥感影像呈现传感器类型多样性、分辨率多样性、时相多样性的特点。这一特点也丰富了海量样本数据成果的内容，有利于智能分类兼顾更多的数据特征。

1.2 高精度地表覆盖分类产品

地理国情监测地表覆盖分类产品包含种植土地、林草覆盖、房屋建筑(区)、铁路与道路、构筑物、人工堆掘地、荒漠与裸露地、水域8个一级类、46个二级类、86个三级类，类型划分以地物自然属性为主要依据[12]，均是基于资源三号、高分二号、WorldView-1/2等高分辨率遥感影像采集，全国范围全覆盖、无缝隙，采集指标基本要求是最小图斑对应的地面实地面积为400 m2。该产品为各部门掌握地表自然资源、生态环境及人类活动基本情况提供了高精度的地表覆盖数据基础。

地理国情监测地表覆盖分类产品格式为矢量空间格式，全国范围图斑数量达到数亿个。在开展样本数据采集时，以县级行政区划为单元，对矢量数据进行了栅格化处理，栅格化粒度为1 m×1 m。采用面积最大法确定各栅格所代表的地表覆盖类型，即栅格范围区域内，各类型中面积占比最大的类型作为该栅格的分类类型。地表覆盖分类产品的分类代码最大为4位数字，因此栅格化结果位深定义为16 bit。

2 研究方法

全国尺度海量样本数据采集数据量巨大，因此，以县级行政区划为单元进行采集，最终聚合形成全国尺度海量样本数据集。基于数据源特征分析，海量样本数据采集需要研究解决县域采集数量权重设置、坐标投影转换、栅格灰度重采样、无效样本数据过滤、地表覆盖分类码映射、样本数据命名标识、特定地表覆盖类型样本数据采集等关键技术问题，实现海量样本数据采集的自动化、批量化，提升采集效率。海量样本数据采集技术流程如图1所示。

2.1 县域采集数量权重设置

每个县级行政区划样本采集数量采用面积因子进行自动分配，根据县域面积与全国县平均面积的大小关系、全国尺度需要采集的样本数据总数量，并考虑采集过程中存在无效样本块的因素，分配合理的采集数量。

2.2 坐标投影转换

地理国情监测地表覆盖栅格数据坐标为CGCS2000，分幅数字正射影像坐标为CGCS2000坐标系下的高斯-克吕格投影坐标。基于位置匹配技术的样本数据采集是由采集的有效的地表覆盖栅格影像块的中心点大地坐标，计算对应的分幅数字正射影像的图幅号，检索分幅数字正射影像数据文件，最后进行采集。

因此，在创建输出的数字正射影像块时，坐标、行列数与地表覆盖栅格影像块保持一致；获取到地表覆盖栅格影像块的中心点大地坐标后，对输出的数字正射影像块进行各像素循环，利用高斯正算算法进行坐标投影转换，从分幅数字正射影像中重采样各像素灰度值。

2.3 栅格灰度重采样

地理国情监测地表覆盖栅格数据分辨率统一为1 m，分幅数字正射影像分辨率包括0.5、1、2 m 3种形式[11]。对于0.5、2 m的分幅数字正射影像，在输出数字正射影像块时，分辨率需要重采样至1 m。

针对分幅数字正射影像存在3种分辨率的情况，在创建输出的数字正射影像块时，坐标、行列数与地表覆盖栅格影像块保持一致；对输出的数字正射影像块进行各像素循环，利用双线性内插算法从分幅数字正射影像中重采样各像素灰度值。

2.4 无效样本数据过滤

2.4.1 数字正射影像无效影像块过滤

采集数字正射影像块时，进行灰度值统计判断，统计像素灰度值为(0,0,0)，即各波段像素灰度值均为0的像素数目。该数目大于影像块像素总数目10%时，将该影像块判定为无效影像块，不进行采集。

在采集过程中，存在另外一种情况，即地表覆盖栅格影像块有效，但数字正射影像块无效。对于有效的地表覆盖栅格影像块，在采集数字正射影像块时，由于数字正射影像按标准分幅组织，可能刚好处于分幅数字正射影像的边缘，因此存在0值区，即数字正射影像块无效的情况。这种地表覆盖栅格影像块不满幅或数字正射影像块不满幅的影像块组合均应删除，以保证样本数据成果的有效性。

2.4.2 地表覆盖栅格无效影像块过滤

采集地表覆盖栅格影像块时，进行灰度值判断，只要存在像素灰度值为0的情况，则说明该影像块存在无值区，不进行采集。

2.5 地表覆盖分类码映射

为便于样本数据应用，并减少数据量，地表覆盖栅格影像块与数字正射影像块均输出为8 bit。由于原始地表覆盖栅格数据为16 bit，因此，需要将其自动转换为8 bit。

转换方法为：创建TXT文本格式的索引文件，索引值为16 bit分类码和对应的8 bit分类码值，针对采集的地表覆盖栅格影像块，根据其原始16 bit分类码，映射得到其8 bit分类码，并作为输出结果值。

2.6 样本数据命名标识

采集的地表覆盖栅格影像块与数字正射影像块均输出为JPG格式，不带空间坐标。数据命名标识规则为：地表覆盖栅格影像块和对应的数字正射影像块文件命名标识一致，由38位字符组成，包含县级行政区划代码、数字正射影像波段数、原始遥感影像获取时间、采集的样本块中心点经度坐标、采集的样本块中心点纬度坐标(用度分秒表示，精度至0.001″)5项信息。具体为：第1—6位字符为6位县级行政区划代码(反映样本数据的空间地域)；第7位字符为数字正射影像波段数；第8—15位字符为原始遥感影像获取时间(从分幅数字正射影像元数据中获取)；第16—25位字符为采集的样本块中心点经度坐标；第26—34位字符为采集的样本块中心点纬度坐标；第35—38位字符为文件后缀名.jpg。例如：1301013201609190980202333360202333.jpg。

从分幅数字正射影像元数据中获取原始遥感影像获取时间时，分为3种情况：①航空摄影影像数据，读取“航摄时间”字段，只有年、月6位字符，第8—15位字符用“0”补齐至8位，如20160900。②卫星影像数据，首选读取“多光谱卫星影像获取时间”字段，有年、月、日8位字符，如20160919；如果该字段值为空，则读取“全色卫星影像获取时间”字段，有年、月、日8位字符，如20160919。③如果读取这3个字段，获取的值均为空，则第8—15位字符用“0”补齐至8位。

2.7 特定地表覆盖类型样本数据采集

样本数据大规模位置匹配采集是以县级行政区划为单元，以地表覆盖栅格数据为主线，通过规则格网划分进行的，那么采集结果就可能出现部分地表覆盖类型采集的样本数目较多、部分地表覆盖类型采集出的样本数目较少甚至没有的情况。同时，部分地表覆盖类型需要采集到纯净的单一地表覆盖类型的样本数据。

针对这种需求，本文采用的方法具体如下：自定义一个TXT文本格式的采集参数文件，将相关信息录入TXT文本文件中，包括批采集数量、特定地表覆盖类型样本数据的分类代码、采集的样本块中心点经度坐标、采集的样本块中心点纬度坐标、样本块列数、样本块行数，各列之间用tab键或空格分开。通过读取采集参数文件，实现相应样本数据的采集，即可满足这种特定需求。

3 软件研发

面向遥感影像智能分类的海量样本数据规模巨大，目前商业软件难以满足采集效率需求。针对这一现状，本文自主研发了栅格空间数据大规模位置匹配采集软件——LCARasterTile，根据地理国情监测高分辨率遥感影像及高精度地表覆盖分类产品的数据管理方式、数据量等特点，按照全国各县面积大小，设置合理的采集参数，利用空间位置匹配技术实现全国尺度海量样本数据采集，提高数据采集效率，满足面向遥感影像智能分类的海量样本数据需要。

3.1 软件特点

LCARasterTile软件是在Visual Studio 2010集成开发环境中，利用Microsoft Visual C++语言实现的。开发过程中，调用了GDAL(geospatial data abstraction library)开源栅格空间数据转换库[13]。

LCARasterTile软件操作简便，实现了大规模栅格空间数据的自动化、批量化位置匹配采集，提高了数据采集效率，保证了样本数据集产品质量，并节约了人力、物力成本。

软件主要具有以下几个特点：①软件支持大数据的输入、读取与处理；②软件可根据遥感影像智能分类对样本影像块尺寸需求，设置不同的采集尺寸；③软件能够实现无效样本数据的自动过滤；④针对特定地表覆盖类型样本数据采集需求，软件可进行特定模式采集；⑤软件具有批处理功能，数据处理效率高；⑥软件界面友好，操作简便。

3.2 软件输入

软件输入数据包括3项：

(1)地理国情监测地表覆盖栅格数据。数据文件采用IMG(后缀名为“.img”)格式，经纬度坐标；多个数据文件存储在同一个文件夹中，便于进行批处理操作。

(2)数字正射影像数据。数据文件采用非压缩的TIFF格式(后缀名为“.tif”)，高斯-克吕格投影坐标；影像坐标信息文件采用TIFF WORLD文档格式(后缀名为“.tfw”)；数据文件与影像坐标信息文件的文件命名标识一致；多个数据文件存储在同一个文件夹中，便于进行批处理操作。

(3)采集参数数据。在进行特定地表覆盖类型样本数据采集时，还需要输入TXT格式采集参数数据；记录内容包括：批采集数量、特定地表覆盖类型样本数据的分类代码、采集的样本块中心点经度坐标、采集的样本块中心点纬度坐标、样本块列数、样本块行数。在进行大规模采集时，采集参数只需输入采集尺寸，默认值为1023，即采集列数、采集行数均为1023个像素。

3.3 软件输出

软件输出成果为位置匹配的遥感影像数据与分类标签数据组成的样本数据对，均为JPG格式，两者命名标识一致，存储在不同的路径下，并一一对应。

4 结果与分析

利用本文方法，基于LCARasterTile软件，以县级行政区划为单元，实现了全国尺度海量样本数据采集。考虑面积、形状等因素，本文选取湖北省武汉市5个县级行政区划作为样例区，利用其成果，对整套方法的实用性及软件运算性能进行评估。

样例区包括武汉市汉南区、蔡甸区、江夏区、黄陂区、新洲区。

采集到的单个位置匹配的遥感影像数据与分类标签数据组成的样本数据对如图2所示，幅面大小为1023 m×1023 m。其中，图2(a)为遥感影像数据(即数字正射影像影像块)，图2(b)为分类标签数据(即地表覆盖栅格影像块)，其像素灰度值代表对应的地表覆盖分类码。

各县级行政区划样本数据采集统计信息见表1，通过采用单个大数据分块、整体数据分批循环处理算法，进行合理的内存分配，解决了数据运行效率与计算性能问题[14]。表1为单台计算机、单线程数据计算效率，在采集过程中，计算机内存占用为300～520 MB(计算机为Windows 7 64位操作系统，32 GB内存，1.9 GHz处理器)。

表1 样本数据采集统计信息

从表1可以看出，县域面积、形状各不同，采集到的样本数据数目、密度存在一定差异。根据表1，采用单台计算机、单线程，全国尺度采集300万量级的样本数据对，需要66 d；如果采用5台计算机、多线程采集，5 d可完成采集，能够满足全国尺度采集效率需求。

5 结 论

(1)本文构建了面向遥感影像智能分类的海量样本数据采集方法、采集软件一整套技术体系，解决了海量样本数据采集中的技术问题。自主研发的LCARasterTile软件，满足了样本数据大规模采集与特定地表覆盖类型采集需求，保证了海量样本数据质量，提升了采集效率。

(2)本文将地理国情监测高分辨率遥感影像及高精度地表覆盖分类产品作为数据源，构建了位置匹配的遥感影像数据与分类标签数据组成的样本数据对，解决了目前高精度且具有标注信息的样本数据不足的问题。

(3)本文利用提出的一整套方法，实现了全国尺度海量样本数据采集，满足了遥感影像智能分类对高质量、大规模样本数据的需求，具有实际意义。采集的全国尺度海量样本数据已应用于遥感影像智能分类研究中，并取得了一定的成果，后续将进一步拓展应用范围。