from
datetime
import
datetime, timedelta
from
satpy.scene
import
Scene
from
pyresample
import
get_area_def
import
warnings
warnings.filterwarnings(
'
ignore
'
)
def
fileNameSplit(fileNs):
fSplit
=
[]
for
fpn
in
fileNs:
fp
=
os.path.dirname(fpn)
fn
= os.path.basename(fpn).split(
'
_P_
'
)[-1
]
nfpn
=
os.path.join(fp, fn)
if
not
os.path.exists(nfpn):
os.rename(fpn, nfpn)
fSplit.append(nfpn)
return
fSplit
def
time_format(info, format=
"
%Y%m%d%H%M%S
"
, format2=
"
%Y%m%d%H%M00
"
):
dt
=
datetime.strptime(info, format)
strdt
= dt + timedelta(hours=8
)
return
datetime.strftime(strdt, format2)
def
get_fy4a(file_paths, channel_list):
fileTime
= time_format(os.path.basename(file_paths[0]).split(
'
_
'
)[-3
])
#
创建scene对象
scn = Scene(file_paths, reader=
'
agri_l1
'
)
#
查看可使用的通道
if
'
C03
'
in
channel_list:
print
(scn.available_dataset_names())
if
'
true_color
'
in
channel_list:
#
查看可用的合成选项
print
(scn.available_composite_names())
#
读取指定通道数据
scn.load(channel_list)
#
打印加载通道后 xarray.Dataset 的数据集
#
print(scn)
>python coord2area_def.py china merc 5 54 60 139 3
### +proj=merc +lat_0=29.5 +lon_0=99.5 +ellps=WGS84
china:
description: china
projection:
proj: merc
ellps: WGS84
lat_0: 29.5
lon_0: 99.5
shape:
height: 2194
width: 2931
area_extent:
lower_left_xy: [-4397119.886334, 553583.846816]
upper_right_xy: [4397119.886334, 7135562.567523]
#
画自定义区域图片
area_id =
'
china
'
area_name
=
"
myarea
"
proj_id
=
'
china_99.5_29.5
'
proj4_args
=
'
+proj=merc +lat_0=29.5 +lon_0=99.5 +ellps=WGS84
'
height
= 2194
width
= 2931
area_extent
= [-4397119.886334, 553583.846816, 4397119.886334, 7135562.567523
]
#
定义地图投影和区域，然后将数据投影到该区域上
areadef =
get_area_def(area_id, area_name, proj_id, proj4_args, width, height, area_extent)
china_scene
=
scn.resample(areadef)
channel_ele
= {
'
C03
'
:
'
vis
'
,
'
C09
'
:
'
vap
'
,
'
C12
'
:
'
ir
'
,
'
true_color
'
:
'
thr
'
}
for
cele
in
channel_list:
pic_dir
= r
'
./
'
pic_name
= f
'
fy4a_{channel_ele[cele]}_{fileTime}.png
'
pic_path
=
os.path.join(pic_dir, pic_name)
if
not
os.path.exists(pic_path):
print
(f
'
生成图片 {pic_name}
'
)
#
保存图片
china_scene.save_dataset(cele, filename=
pic_path)
if
__name__
==
'
__main__
'
:
t1
=
datetime.now()
grayFilePathList
= [
'
./FY4A-_AGRI--_N_REGC_1047E_L1-_FDI-_MULT_NOM_20211020022336_20211020022753_4000M_V0001.HDF
'
,
'
./FY4A-_AGRI--_N_REGC_1047E_L1-_GEO-_MULT_NOM_20211020022336_20211020022753_4000M_V0001.HDF
'
]
print
(grayFilePathList)
channelList
= [
'
C03
'
,
'
C09
'
,
'
C12
'
]
get_fy4a(grayFilePathList[:
1
], channelList)
t2
=
datetime.now()
print
(t2 -
t1)
print
(
'
------
'
)
thrFilePathList
= [
'
./FY4A-_AGRI--_N_REGC_1047E_L1-_FDI-_MULT_NOM_20211020022336_20211020022753_4000M_V0001.HDF
'
,
'
./FY4A-_AGRI--_N_REGC_1047E_L1-_GEO-_MULT_NOM_20211020022336_20211020022753_4000M_V0001.HDF
'
]
print
(thrFilePathList)
channelList1
= [
'
true_color
'
]
get_fy4a(thrFilePathList, channelList1)
t3
=
datetime.now()
print
(t3 - t2)
reader='agri_l1'
版本<=0.36.0 'agri_l1' >=0.37.0 'agri_fy4a_l1' 其他请看网址
https://satpy.readthedocs.io/en/stable/index.html#reader-table
HTTPSConnectionPool(host='zenodo.org', port=443): Max retries exceeded with url: /record/1288441/files/pyspectral_atm_correction_luts_no_aerosol.tgz
# 原因下载依赖文件失败
# 手动下载地址
https://zenodo.org/record/3824535/files/pyspectral_rsr_data.tgz
https://zenodo.org/record/1288441/files/pyspectral_atm_correction_luts_no_aerosol.tgz
pyspectral_rsr_data.tgz 解压放到下面地址
C:\Users\admin\AppData\Local\pytroll\pyspectral
pyspectral_atm_correction_luts_no_aerosol.tgz 解压放到下面地址
C:\Users\admin\AppData\Local\pytroll\pyspectral\rayleigh_only
linux
pyspectral_rsr_data.tgz 解压放到下面地址
/root/.local/share/
pyspectral
pyspectral_atm_correction_luts_no_aerosol.tgz 解压放到下面地址
/root/.local/share/pyspectral/rayleigh_only
Don't know how to open the following files: {'./FY4A.......HDF'}
# 按我猜测，应该是satpy 和其依赖包版本的问题
h5py==3.7.0
numpy==1.21.6
satpy==0.36.0
pyresample==1.22.0
scikit-image==0.19.3
pyspectral==0.12.0
