张 镇,李水平,张志远
(1.武汉大学 遥感信息工程学院,湖北 武汉 430079;2.合肥工业大学 土木与水利工程学院,安徽 合肥 230009;3.中国科学院精密测量科学与技术创新研究院 大地测量与地球动力学国家重点实验室,湖北 武汉 430077;4.水利部信息中心,北京 100053)
遥感技术是现代信息技术发展的重要成果之一,是一门新兴实用的、先进的、不断发展的空间探测技术[1-3];而遥感影像是用户与遥感技术研究、应用之间的桥梁,是获取地物信息的基础数据,具有实时性高、覆盖范围广、信息丰富客观等优点[4]。卫星遥感技术的发展使得遥感影像数据量迅速增加,呈现出全球、海量、多源的特征[5-8]。查询满足特定时间范围、光谱波段、空间分辨率要求的遥感数据,并保证多源遥感影像的覆盖率满足应用需求,是大范围遥感业务应用首先需要面对的一个问题[9]。抗洪、抗震等应急救灾工作更需要快速查询满足覆盖率要求的高分辨率卫星遥感影像,进而利用这些数据获取灾区现状和评估信息[10]。传统的卫星影像覆盖率获取方法主要是先在相关平台手动查询、搜索和对比对应的卫星遥感影像,再通过影像快视图查看整体的覆盖情况。整个过程需在平台界面进行交互,影像及其覆盖率查询过程较复杂、效率较低,不利于项目和业务的进展。苗立新[9]等利用ENVI/IDL实现了多源遥感影像覆盖范围快速查询技术,可查询总体的覆盖情况,但不能获取精确的影像覆盖率值;我国自然资源卫星遥感云服务平台提供了各省级行政区划每年的实时遥感影像覆盖率统计模块,但用户无法获取目标区域某个较小时间段内的影像覆盖率信息,因此需要一种可根据需求快速获取精确的卫星遥感影像覆盖率值的方法,以满足遥感技术日益发展的需求。
为方便用户共享遥感影像资源,相关机构推出了开源遥感影像数据访问API,允许用户通过编程语言访问平台中存储的数据,进而获取所需的卫星遥感影像。基于这一技术背景,本文以开源遥感数据访问API为技术支撑,利用Python相关模块编程实现了卫星遥感影像覆盖率查询方法,能更加方便地获取卫星影像覆盖率信息,并以目标区域覆盖率实时查询监测为应用案例进行了分析和总结。
在遥感技术中,遥感影像的应用水平明显滞后于航天遥感技术的发展,使得很大一部分卫星影像资源处于搁置状态,形成了数据生产以及传输能力远超数据获取、处理、应用能力的局面[11]。为方便用户共享卫星遥感数据,相关机构推出了开源遥感影像数据访问API,允许用户通过编程语言查询卫星影像数据,使得获取遥感数据更加方便,从而推动遥感影像数据的应用。本文对国内外相关开源遥感数据访问API及其主要信息进行了总结,可以更清楚地了解这些应用程序接口的内容和特点,如表1所示。
表1 开源遥感数据访问API
本文以上述遥感数据访问API为基础,以Python为编程开发语言进行相关处理与分析,实现了卫星遥感影像覆盖率的快速查询。技术流程如图1所示,首先确定必要的查询参数,并利用上述API查询所需卫星影像的元数据;再通过Python相关模块对影像元数据进行一系列的转换、地理处理与分析;最后计算得到准确的覆盖率值。
图1 卫星遥感影像覆盖率查询流程图
首先根据开源遥感数据访问API文档中的内容和功能构造特定查询参数,再利用相关工具发送查询请求,并获取服务器响应数据,即可查询所需卫星遥感影像的元数据。本文利用Python相关模块设计程序实现了卫星遥感影像元数据的查询功能,具体步骤为:
1)构建查询参数。与在平台查询影像相同,利用遥感数据访问API查询遥感数据也需要确定查询参数。遥感数据访问API规定了空间范围、时间范围、卫星和传感器种类等必要参数,EarthExplorer API和OpenSearch API还包括月份限定、查询数量限制、排序方法等其他参数,如表2所示。本文利用Python相关方法,将这些参数构建成标准的、特定数据格式的查询参数变量,作为发送查询请求的输入参数。
表2 查询参数汇总
2)发送查询请求。将构建的查询参数作为输入参数,向服务器发送查询请求并获取响应数据,即可获得遥感影像的元数据。本文利用Python中的urllib2和requests模块发送查询请求。urllib2是Python中的一个 HTTP客户端库[12-13],为内置模块,其request()方法可构造一个请求对象,将其作为该模块urlopen()方法的输入参数即可向服务器发送查询请求。requests也是一个HTTP库[14],为第三方模块,需单独安装,其在Python内置模块的基础上进行了高度封装,可实现网页和浏览器请求的所有操作。url为查询请求网址,jsonData为查询参数变量,利用urllib2和requests模块发送查询请求的核心代码为:
request = urllib2.Request(url=url, data=jsonData, headers=header) ‘ urllib2发送查询请求
response = urllib2.urlopen(request)
response = requests.get(url, params=jsonData, headers=header) ‘ requests发送查询请求
3)存储与展示数据。遥感影像元数据是元数据的一种,通过元数据建立海量卫星遥感数据目录和数据交换中心,为空间数据的共享提供了可能[15]。以通过资源卫星中心API查询得到的影像元数据为例,其数据格式为JSON[16-17],本文首先对其进行解析处理;然后利用Python中的JSON模块将其存储为后缀名为json的文件[18],并通过JSON在线解析工具(https://www.json.cn/)浏览查询结果。元数据内容主要包括卫星影像的角点坐标、对应影像的空间覆盖范围以及行列号、卫星ID、成像时间、传感器、云覆盖量、产品号等其他元数据信息。经过解析处理的一景GF-2号卫星影像元数据示例为:
’影像元数据示例
{ ’其他元数据信息
“PathRow”: “1010165”,
“satellite”: “GF2”,
“title”: “GF2_PMS1_4088967”,
“imgRectangle”: Array[5],
“sceneId”: 6594578,
“date”: “2019/07/02”,
“icon”:“http://218.247.138.119:7777/DSSPlatform/BROWSE/ARCHIVE_PGS/GF2/LEVEL1A/PMS1/20190702/GF2_PMS1_E117.1_N30.3_20190702_L1A0004088967.jpg”,
“sensor”: “PMS1”,
“cloud”: 0,
“productId”: 4088967
}
“imgRectangle”: [ ’角点坐标信息
[117.037, 30.418],
[117.273, 30.3734],
[117.221, 30.1693],
[116.986, 30.2139],
[117.037, 30.418]
]
本文利用Python相关功能模块对查询获取的卫星影像元数据进行了一系列处理与分析,最终计算得到卫星影像覆盖率值,进而实现了卫星影像覆盖率查询功能。
1)将影像元数据转换为GeoJSON。利用Python中的GeoJSON模块将每景影像的元数据转换为对应的面状GeoJSON数据。GeoJSON是一种对各种地理数据结构进行编码的格式,支持点、线、面等几何类型,其结构包含一个或多个要素,每个要素又包含一个 几何对象和其他属性两个部分[19-20]。提取元数据中的角点坐标作为坐标数据构成GeoJSON几何对象中的坐标部分,其他元数据信息则构成GeoJSON的属性信息部分。具体转换流程如图2所示。
图2 影像元数据转换为GeoJSON的流程图
利用§2.1中的影像元数据示例转换得到的GeoJSON数据为:
“properties”:{ ‘ 属性信息
“sceneId”: 6594578,
“satellite”: “GF2”,
“title”: “GF2_PMS1_4088967”,
“imgRectangle”: Array[5],
“PathRow”: “1010165”,
“date”: “2019/07/02”,
“productId”: 4088967,
“sensor”: “PMS1”,
“cloud”: 0,
“icon”:“http://218.247.138.119:7777/DSSPlatform/BROWSE/ARCHIVE_PGS/GF2/LEVEL1A/PMS1/20190702/GF2_PMS1_E117.1_N30.3_20190702_L1A0004088967.jpg”
}
“geometry”:{ ‘ 几何对象
“type”: “Polygon”,
“coordinates”: [
[
[117.037, 30.418],
[117.273, 30.3734],
[117.221, 30.1693],
[116.986, 30.2139],
[117.037, 30.418]
]
]
}
2)将GeoJSON转换为面状矢量要素。利用Python中的GeoPandas模块将所有面状GeoJSON转换为具有指定地理坐标的面状矢量要素。其中,GeoJSON中的几何对象对应的是矢量要素的空间几何属性,即遥感影像的空间范围;GeoJSON中的属性信息则存储在面状矢量的属性表中。具体转换流程如图3所示,GeoJSON转换的面状矢量要素属性表中的部分信息如图4所示。该方法在转换时直接将元数据信息添加到了要素属性表中,因此可在ArcGIS中同时查看影像覆盖情况和影像元数据信息。
图3 GeoJSON转换为面状矢量要素流程图
图4 面状矢量要素属性表部分信息
3)计算覆盖率。ArcPy是一套ArcGIS处理脚本的软件包,其目的是利用Python建立实用且高效的地理数据分析、转换、管理以及自动制图的处理工 具[21-22]。Python中的ArcPy模块包括投影、联合、相交、计算面积等与ArcGIS中相关功能对应的方法,具体代码为:
arcpy.Project_management() ‘投影
arcpy.Union_analysis() ‘联合
arcpy.Intersect_analysis() ‘相交
arcpy.CalculateAreas_stats() ‘计算面积
首先利用arcpy.Project_management()方法对所有面状矢量进行批量投影;再利用arcpy.Union_analysis()方法将投影后的面状矢量联合为一个图层,这样可避免对卫星影像重叠区域面积的重复计算,保证计算结果的正确性;然后利用arcpy.Intersect_analysis()方法将联合后的整体矢量要素与具有相同投影坐标的目标区域矢量图层相交;最后利用arcpy.CalculateAreas_stats()方法计算相交结果的面积。具体操作流程如图5所示。
图5 覆盖率计算部分技术流程图
卫星遥感影像覆盖率为相交图层的面积与目标区域矢量图层面积的比值,计算公式为:
式中,C为覆盖率值;S1为相交结果图层的面积;S2为目标区域矢量图层的面积。
利用上述卫星遥感影像覆盖率查询方法,只需明确必要的查询参数,即可快速获取所需遥感影像数据及其对目标区域的覆盖率值,具有实时、准确度高等优势。若覆盖率不满足应用要求,则需寻找其他卫星影像数据源作为补充。计算覆盖率时,由于需要对大量矢量数据进行批量处理,当影像数量较多时,程序运行较慢,需进一步改进算法,以提高程序运行效率。
我国是一个地震、泥石流、洪涝等自然灾害频发的国家,重大自然灾害大多具有突发性强、涉及面广、危害性大等特点,因此应急救灾工作需要利用满足覆盖率要求的高分辨率卫星影像获取灾区现状和评估信息[10]。高分辨率卫星影像的覆盖率能否满足要求将直接影响应急救灾工作的进展,而目标区域的卫星影像覆盖率也反映了该地区卫星影像的获取和应用能力。因此,实时查询和掌握我国各地区的高分辨率卫星遥感影像覆盖率情况,对于灾害应急响应领域具有重要意义。
夏季是泥石流、洪涝等灾害的多发期,灾害一旦发生,需及时查询和分析灾害发生前后的遥感影像和其他相关数据,从而掌握各种设备、设施以及人民群众可能的受损情况,为抗灾救援决策提供及时的数据支持。鉴于此,本文以2019年7月我国各省级行政区划1 m和2 m分辨率的卫星影像为处理对象,进行了覆盖率分区域查询和监测分析;并以此为应用案例,开展了实验与分析工作。
由卫星遥感影像覆盖率查询流程可知,完成这样的实验需要根据需求确定必要的查询参数作为输入数据,才能查询卫星影像覆盖率值。本文设计的查询方案如图6所示。
图6 覆盖率查询方案
本文对查询得到的各省级行政区划(由于澳门特别行政区面积太小,查询意义不大,故省略)1 m和 2 m分辨率的卫星影像覆盖率值进行汇总,结果如表3所示(某些省级行政区划名称使用了简称)。
表3 2019年7月全国省级行政区划1 m和2 m分辨率 卫星影像覆盖率统计表
为总体展示全国各地卫星遥感影像的覆盖情况,本文利用ArcGIS将覆盖率作为属性值添加至全国省级行政区划矢量要素的属性表中,并分别创建1 m和 2 m分辨率卫星遥感影像覆盖率柱状图,如图7、8所示,图中颜色渐变表示覆盖率的高低。
卫星遥感影像覆盖率的高低反映了该地区遥感影像的获取和应用能力,由表3和图7、8可知,在特定时间段、特定区域内,2 m分辨率卫星遥感影像覆盖率普遍比1 m分辨率卫星遥感影像高;且在较短时间内,1 m分辨率卫星影像数量较少,覆盖率较低,有的地区甚至没有影像数据,存在覆盖率为零的情况,因此需要进一步提升高分辨率卫星影像的生产能力以及相关地区的影像获取和应用能力。从总体来看,无论是1 m分辨率还是2 m分辨率卫星遥感影像,均在我国中东部等高纬度地区分布数量较多,覆盖率值较高,其他地区数量较少、覆盖率值较低,由此可知我国高分辨率卫星影像的数量在空间上分布并不均匀,需在提升卫星影像生产能力的同时做到全国区域全覆盖。
图7 2019年7月省级行政区划1 m分辨率卫星遥感影像覆盖率柱状图
图8 2019年7月省级行政区划2 m分辨率卫星遥感影像覆盖率柱状图
由上述应用案例分析可知,通过本文方法可实时查询和监测目标区域特定分辨率卫星影像的覆盖率信息,进而了解该地区卫星遥感影像的生产、获取、应用情况,可为灾害应急、国土规划等工作提供一定的数据支持。
针对传统卫星遥感影像覆盖率信息获取途径较繁琐、获取时间较长、无法获取准确值等问题,本文提出了基于遥感数据访问API的卫星遥感影像覆盖率查询方法。该方法具有的优势为:
1)仅需确定必要的查询参数,即可查询获取对应的卫星遥感影像数据及其对目标区域的覆盖率信息,比传统方法更高效、准确,极大地节省了时间。
2)既可根据实际应用需求查询所需的卫星影像数据及其覆盖率值,又可实时查询和监测目标区域特定分辨率卫星影像的覆盖率,从而了解目标区域卫星遥感数据的获取和应用情况。
本文方法降低了获取遥感影像数据及其覆盖率信息的时间成本,可为某些遥感相关技术领域提供数据和技术支持,具有一定的应用价值;但该方法计算覆盖率时需对大量矢量数据进行批量处理,读取、输出数据将消耗较多时间,降低了程序运行效率,因此需进一步优化算法,以提高程序运行效率。