基于python包sentinelsat的Sentinel-2数据下载

关于Sentinel-2数据下载和预处理的介绍,网上已经有很多很好的资源,可以参考超级禾欠水 的系列文章。

sentinelsat包虽然提供了下载数据的功能,但是有时候受网速的影响,下载可能会中断,这时候重新运行代码,下载会重新开始,暂时还没有找到断点续传的功能。因此,这里介绍一种先利用sentinelsat包查询Sentinel-2图像,基于查询信息生成 wget 命令,然后利用 wget 下载Sentinel-2图像的方法。

sentinelsat 查询数据信息

获得研究区域的遥感影像,通常有两种方法,一种是知道研究区域的经纬度范围,在查询系统中选择区域查询数据,另一种是知道覆盖研究区影像的编号,例如Landsat的path和row,Sentinel-2 的MGRS tile编号,通过编号查询影像。sentinelsat也提供了两种查询方式。

直接上代码

from collections import OrderedDict
from sentinelsat import SentinelAPI, read_geojson, geojson_to_wkt
import requests

def query_geo(api, footprint, date, producttype, cloudcover):
    # 通过设置OpenSearch API查询参数筛选符合条件的所有Sentinel-2 L2A级数据
    products =api.query(footprint,        # Area范围
            date=date, # 搜索的日期范围
            platformname='Sentinel-2',    # 卫星平台名,Sentinel-2
            producttype=producttype,        # 产品数据等级,‘S2MSI2A’表示S2-L2A级产品
            cloudcoverpercentage=cloudcover)  # 云量百分比
    return products

def query_MGRS(api, tiles, date, producttype, cloudcover):
    # tiles = ['50SKE']
    query_kwargs = {
        'platformname': 'Sentinel-2',
        'producttype': producttype,
        'date': date}

    products = OrderedDict()
    for tile in tiles:
        kw = query_kwargs.copy()
        kw['filename'] = '*_T{}_*'.format(tile)  # products after 2017-03-31
        pp = api.query(cloudcoverpercentage=cloudcover, **kw) # 设置云量筛选是有效的
        products.update(pp)

    return products

if __name__ == "__main__":
    # 请使用哥白尼数据开放获取中心自己的用户名及密码
    username = '*****'
    password = '*****'

    cmd_str = r'wget -c --http-user={} --http-password={} --no-check-certificate -O '.format(username, password)
    out_dir = r'E:\temp' + '\\'
    cmd_str += out_dir
    bat = open(out_dir + 'wget_download.bat', 'w')

    date = ('NOW-14DAYS', 'NOW')
    # date = ('20190101', '20191231')
    # producttype: Sentinel-2: S2MSI2A, S2MSI1C, S2MS2Ap
    producttype = 'S2MSI1C'
    # cloudcoverpercentage 云覆盖筛选
    cloudcover = (0, 75)
    # 创建SentinelAPI
    api = SentinelAPI(username, password, 'https://scihub.copernicus.eu/dhus//')

    # 两种不同的查询方式
    # -------------------------------------------------------------
    # 1. 利用footprint,两种得到footprint的方法
    # (1)直接定义footprint
    # footprint = 'POLYGON((114.1631 35.5892,114.5792 35.5892,114.5792 35.8746,114.1631 35.8746,114.1631 35.5892))'
    # (2)绘制geojson,读取得到footprint,在 http://geojson.io 绘制保存为 geojson
    # geojson = r'E:\temp\temp.geojson'
    # 读入geojson文件并转换为wkt格式的文件对象
    # footprint = geojson_to_wkt(read_geojson(geojson))

    # products = query_geo(api, footprint, date, producttype, cloudcover)
    # -------------------------------------------------------------
    # 2. 利用MGRS编号
    tiles = ['50SKE', '49SGT']
    products = query_MGRS(api, tiles, date, producttype, cloudcover)
    # -------------------------------------------------------------

    #通过for循环遍历查询到的产品文件名
    for product in products:
        # 通过OData API获取单一产品数据的主要元数据信息
        # 包含产品的id, title, size, url, md5sum(校验码), date
        product_info = api.get_product_odata(product)
        # print(product_info.keys())
        # 打印下载的产品数据文件名
        fname = product_info['title']
        print(fname)
        url = product_info['url']
        # print(product_info['date'])
        # print('id: {}'.format(product_info['id']))
        # print('md5: {}'.format(product_info['md5']))

        # # 下载数据的快视图
        jpgname = fname + '-ql.jpg'
        qlook_url = url[0:(-len('$value'))] + "Products('Quicklook')/$value"
        # print(qlook_url)
        # 需要设置登录信息,否则status_code为401
        # 401错误表示的是被请求的页面需要用户名和密码
        response = requests.get(qlook_url, stream=True, auth=(username, password))
        if response.status_code == 200:
            with open(out_dir + jpgname, 'wb') as f:
                # 每128个流遍历一次
                for data in response.iter_content(128):
                    # 把流写入到文件
                    f.write(data) # data相当于一块一块数据写入到我们的图片文件中
        else:
            print(jpgname + ' 快视图下载失败')

        # 生成wget下载的命令
        bat.write(cmd_str + fname + '.zip ' + url +'\n')
        # 下载产品id为product的产品数据
        # api.download(product)

    bat.close()

代码中query_geo是基于坐标范围查询的函数,其中footprint是区域的经纬度坐标,获得经纬度,一方面可以可以在geojson里面绘制并保存研究区的geojson格式文件,在代码中通过footprint = geojson_to_wkt(read_geojson(geojson))读取;另一方面,还可以Sentinel数据下载网站里, 查询数据后,从设置的查询条件中得到,如下图:
在这里插入图片描述对应的查询条件为
Request Done: ( footprint:“Intersects(POLYGON((114.11367826359871 35.839385438727206,114.23972512137068 35.839385438727206,114.23972512137068 35.93066306812729,114.11367826359871 35.93066306812729,114.11367826359871 35.839385438727206)))” ) AND ( beginPosition:[2020-03-01T00:00:00.000Z TO 2020-03-14T23:59:59.999Z] AND endPosition:[2020-03-01T00:00:00.000Z TO 2020-03-14T23:59:59.999Z] ) AND ( (platformname:Sentinel-2 AND producttype:S2MSI1C))

query_MGRS可以通过MGRS查询数据,输入的tiles是需要查询影像的MGRS列表,MGRS可以从文件名中得到:
S2A_MSIL1C_20200303T030621_N0209_R075_T49SGT_20200303T061730

代码运行结束后,把wget命令保存在文件 wget_download.bat中,可以用记事本打开,并可以在cmd命令行中运行。

代码还有一个功能,就是可以下载快视图,可以先通过查询研究区数据,并下载快视图,通过快视图,筛选出需要的影像,然后再修改保存的wget命令文件 wget_download.bat,进行文件下载

把下载的 wget.exe 放在 out_dir 文件夹下面,在命令行中运行wget_download.bat开始执行下载。
在这里插入图片描述

  • 4
    点赞
  • 12
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值