支持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/
有以下三种:
原始数字量化值 (所有通道)
反射率 (C01 - C06)
辐射率和亮温 (C07 - C14)
$ 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。