风云4A(FY4A) AGRI L1数据读取基础介绍¶

1. 支持的文件类型¶

支持NetCDF格式的风云4A(FY4A) AGRI数据(全圆盘和中国区)

示例:

全圆盘:

FY4A-_AGRI--_N_DISK_1047E_L1-_FDI-_MULT_NOM_20190807060000_20190807061459_4000M_V0001.HDF

中国区:

FY4A-_AGRI--_N_REGC_1047E_L1-_FDI-_MULT_NOM_20190807045334_20190807045750_1000M_V0001.HDF

全圆盘标识：DISK , 中国区标识：REGC.

数据链接:

实时数据(30天内)及文档:
    https://fy4.nsmc.org.cn/data/en/data/realtime.html
历史数据 (2018-03-12 -- ):
    http://satellite.nsmc.org.cn/PortalSite/Data/Satellite.aspx
FY4A官方应用平台:
    http://rsapp.nsmc.org.cn/geofy/

2. 定标¶

有以下三种:

原始数字量化值 (所有通道)
反射率 (C01 - C06)
辐射率和亮温 (C07 - C14)

3. 例子¶

安装¶

$ conda install -c conda-forge satpy

加载数据¶

import os, glob
from satpy.scene import Scene

# 加载FY4A文件
filenames = glob.glob('/xin/data/FY4A/20190807/FY4A-_AGRI*4000M_V0001.HDF')

# 创建scene对象
scn = Scene(filenames, reader='agri_l1')

# 查看可用的通道
scn.available_dataset_names()

['C01',
 'C02',
 'C03',
 'C04',
 'C05',
 'C06',
 'C07',
 'C08',
 'C09',
 'C10',
 'C11',
 'C12',
 'C13',
 'C14',
 'satellite_azimuth_angle',
 'satellite_zenith_angle',
 'solar_azimuth_angle',
 'solar_glint_angle',
 'solar_zenith_angle']

# 以红外通道为例
ir_channel = 'C12'
scn.load([ir_channel], generate=False, calibration='brightness_temperature')

全圆盘图¶

# 在notebook中显示
scn.show(ir_channel)

# 保存到文件
# scn.save_dataset(ir_channel, filename='{sensor}_{name}.png')

真彩色全圆盘图¶

# 查看可用的合成选项
scn.available_composite_names()

Too many possible datasets to load for DatasetID(name=None, wavelength=3.9, resolution=None, polarization=None, calibration=None, level=None, modifiers=())
Too many possible datasets to load for 3.9
Too many possible datasets to load for 3.9
Too many possible datasets to load for 3.9
Too many possible datasets to load for DatasetID(name=None, wavelength=3.9, resolution=None, polarization=None, calibration=None, level=None, modifiers=())
Too many possible datasets to load for 3.9
Too many possible datasets to load for 3.9

['ash',
 'dust',
 'fog',
 'green',
 'green_snow',
 'ir108_3d',
 'ir_cloud_day',
 'natural_color',
 'natural_color_sun',
 'night_background',
 'night_background_hires',
 'overview',
 'overview_sun',
 'true_color']

# 注：这步需要大内存 (取决于cpu核数)
# 可以查看FAQ中关于内存的讨论：
#    https://satpy.readthedocs.io/en/latest/faq.html

composite = 'true_color'
scn.load([composite])
scn.show(composite)
# scn.save_dataset(composite, filename='{sensor}_{name}.png')

Required file type 'agri_l1_4000m_geo' not found or loaded for 'solar_zenith_angle'
Required file type 'agri_l1_4000m_geo' not found or loaded for 'satellite_azimuth_angle'
Required file type 'agri_l1_4000m_geo' not found or loaded for 'satellite_zenith_angle'
Required file type 'agri_l1_4000m_geo' not found or loaded for 'solar_azimuth_angle'

特定区域图¶

以下以台风利奇马为例。

首先，我们需要定义地图投影和区域，然后将数据投影到该区域上。

我们可以用Pyresample来定义区域，并将其写入到area.yaml方便以后直接调用。

from pyresample import get_area_def

area_id = 'lekima'

x_size = 549
y_size = 499
area_extent = (-1098006.560556, -967317.140452, 1098006.560556, 1026777.426728)
projection = '+proj=laea +lat_0=19.0 +lon_0=128.0 +ellps=WGS84'
description = "Typhoon Lekima"
proj_id = 'laea_128.0_19.0'

areadef = get_area_def(area_id, description, proj_id, projection,x_size, y_size, area_extent)

你可以用coord2area_def.py程序来直接生成区域定义。

比如用python coord2area_def.py lekima_4km laea 10 28 118 138 4可得到之前定义的利奇马区域：

lekima_4km:
  description: lekima_4km
  projection:
    proj: laea
    ellps: WGS84
    lat_0: 19.0
    lon_0: 128.0
  shape:
    height: 499
    width: 549
  area_extent:
    lower_left_xy: [-1098006.560556, -967317.140452]
    upper_right_xy: [1098006.560556, 1026777.426728]

然后将该定义拷贝到$PPP_CONFIG_DIR/areas.yaml中，即可直接调用。

# 如果你已经添加区域到areas.yaml中，可直接调用:
os.environ['PPP_CONFIG_DIR'] = '/yin_raid/xin/satpy_config/'
lekima_scene = scn.resample('lekima_4km')

# 否则需要使用之前定义的区域:
# lekima_scene = scn.resample(areadef)

lekima_scene.show(composite)
# lekima_scene.save_dataset(composite, filename='{sensor}_{name}_resampled.png')

如果想利用自定义的colormap来生成图像（如下图），请参阅关于enhancement的notebook。