基于官方API实现ESA遥感数据批量下载的方法

本文链接：https://blog.csdn.net/Jcb1906824038/article/details/125984139

一、确定Sentinel-2影像瓦片Tile分区号

Sentinel-2 A/B卫星的观测范围在56°S-84°N之间，影像数据采用UTM系统分区编码规则。其中，每个100km×100km的瓦片（Tile）完整的编号需要由”带区“标识和网格”行列“组合确定。

“带区标识”：将地球表面按纬差8°、经差6°分成”带区“，从经度180°往东，每个经差6°的区域以1-60的带号进行编号；从80°S往北算起，每8°纬差分别按照C-X（去掉I、O）20个字母顺序进行区号表示（具体区号标识图解如下所示）。

“网格行列”：在”带区“网络区域标识的基础上，进一步划分出边长为100km的网格单元。网格单元编号由分别代表东西方向列号和南北方向行号的2个英文字母组合表示。

东西方向上，列号从1带的西边缘开始用A-Z（去掉I、O）24个字母依次反复表示，每3个带字母反复一次，每带8个字母，对应完整的或不完整的8个窄带（在赤道附近，每带共有6个完整和2个不完整的窄带）。

南北方向上，奇数投影带的行号从赤道开始用A-V（去掉I、O）18个字母往北依次反复表示，约经过2个多带区字母反复一次。偶数投影带的行号从赤道以南的第5个网格处开始用A-V（去掉I、O）18个字母往北依次反复表示，这样赤道以北第1个100km网格就是F。这样奇数和偶数带相差5个字母可以错开相同代号，避免指示目标和读图时发生混乱。

举个例子：（确定覆盖武汉的Sentinel-2遥感影像的Tile编号的方法）

下载武汉市的行政区划矢量文件并在ESA官网获取Sentinel-2 1C数据的kml格式平铺网格(详见参考资料[10])。首先，利用Arcgis系统转换工具将kml格式文件转换为图层格式；然后，一般情况下，将行政区划矢量边界地理坐标系转为GCS_WGS_1984与平铺网格一致（或根据具体项目要求转换合理坐标系）。最后，两者叠加观察并打开平铺网格的属性表观察name字段即可。

二、Python利用sentinelsat(结合urllib3库)、collections库包进行自动获取

（一）首次在Pycharm中使用Python运行sentinelsat库代码会出现virtualenv无法加载文件activate.ps1的情况，这是因为Windows10系统默认执行的是严格受限模式。解决方式如下：(若需进一步学习，可详见参考资料[15][16])

以管理员身份打开PowerShell执行Set-ExecutionPolicy Bypass命令，然后重启pycharm即可。

（二）以一种奇怪的方式解决了SSL版本问题

在运行代码时，需要关闭科学上网，否则会出现SSL问题。这个问题我认为可能出自代理，后续值得研究一下。

（二）代码块

此处选择采用https://scihub.copernicus.eu/dhus而不是默认值，是因为默认值有可能无法获取早期遥感影像数据[8]。

from sentinelsat import SentinelAPI
from collections import OrderedDict
import urllib.request
import urllib3

# 检查网站接入情况，网站状态码为200表示正常
response = urllib.request.urlopen('https://scihub.copernicus.eu/dhus')
if response.status != 200:
    print(response.status)
    print("网站接入异常")
    pass
else:
    print("网站成功接入")
api = SentinelAPI('username', 'password', 'https://scihub.copernicus.eu/dhus')
# api接入
if api:
    print(api)
else:
    print("api接入失败")
# 此处将需要下载的瓦片编号组成列表存储
tiles = ['49RGP', '49RGQ', '50RKU', '50RKV']
query_kwargs = {
    'platformname': 'Sentinel-2',
    'producttype': 'S2MSI1C',
    'date': ('20200520', '20200620'),
    'filename': 'S2A_*'
               }
products = OrderedDict()
for tile in tiles:
    kw = query_kwargs.copy()
    kw['tileid'] = tile
    pp = api.query(**kw)
    products.update(pp)
    print(tile)
# 显示满足此要求下的影像数量
print(len(products))
api.download_all(products)


# 另一种根据roi选择影像的方法
'''
api = SentinelAPI('username', 'password', 'http://scihub.copernicus.eu/dhus')
roi = 'POLYGON(( 38.34503173828125 29.826646200964564,\
                  38.356361389160156 29.826646200964564,\
                  38.356361389160156 29.836325670014148,\
                  38.34503173828125 29.836325670014148,\
                  38.34503173828125 29.826646200964564))'
products = api.query(area=roi, date=('20190524', '20190530'), producttype='S2MSI2A', cloudcoverpercentage=[0, 50])
downfiles = OrderedDict()
for i in products:
    product = products[i]
    filename = product['filename']
    print(filename)
'''