☰
FY4A数据读取代码包:从HDF5到定标投影的完整链路
2026/10/16 1:43:16 网站建设 项目流程

简介:这份资源是面向气象遥感、卫星数据处理方向的开发者与科研人员整理的FY4A数据读取代码包,针对风云四号静止气象卫星高分辨率观测数据的解析与二次开发需求,提供了一套可直接运行的Python实现框架。压缩包共6个文件,以5个py脚本和1个xml配置为主,整体约12KB,体量轻巧但结构完整:主程序负责调度读取流程,工具模块承担数据转换与预处理,组件层封装解码与图像处理逻辑,配置项则集中管理数据源路径与输出参数,另附测试脚本用于验证功能正确性。目前已有2299人学习下载,说明其在同类资料中具备一定参考价值。读者可借此快速理解FY4A原始二进制数据的结构化解析思路,掌握从数据读取、校验到分析输出的完整链路,并在此基础上搭建自己的气象监测或环境灾害预警处理系统,适合作为入门学习与项目原型开发的参考脚本。

1. 拆开这份 FY4A 数据读取代码包:它到底能帮你省掉哪几天

如果你手头有一批风云四号 A 星(FY4A)的 HDF5 级数据,却卡在“怎么把文件读进内存、怎么把定标和投影做对”这一步,这份FY4A数据读取代码codes.7z就是冲这个场景来的。它不是一个完整业务系统,而是一套围绕 FY4A 数据读取、通道拆分、经纬度定位和基础可视化的小型代码集合,适合做遥感数据处理、气象算法验证、卫星产品预研的从业者直接拿去改。风云四号代码这类资源网上散落不少,但多数只给一段读文件的片段,缺定标、缺投影、缺批量处理,拿过来跑不通。这份包的价值在于把“读—定标—定位—出图”这条链路串起来了,新手能照着跑通第一张图,熟手能直接抠出 AGRI 定标表和投影参数去接自己的流程。下面我按实际拆包顺序讲,先说什么文件对应什么功能,再讲怎么改参数、怎么避坑。

2. FY4A 数据格式与代码包结构:先搞清 HDF5 里装了什么

2.1 FY4A 二级产品的基本组织方式

FY4A 的 AGRI 载荷数据一般以 HDF5 格式分发,常见的是 L1 级全圆盘或区域产品。一个文件里通常包含多个数据集(dataset),核心是NOMChannel01到NOMChannel14这类通道数据,外加CALChannel01到CALChannel14定标表、NOMSatelliteZenith等角度信息,以及文件级的全局属性(global attributes)记录投影参数和观测时间。很多人第一次用h5py打开文件,看到一堆 key 就懵了,不知道哪个是辐亮度、哪个是亮温、哪个是反射率。实际上 FY4A 的可见光通道(1-6)和红外通道(7-14)定标方式不同:可见光走反射率定标,红外走亮温定标,定标表里存的是查找表(LUT),不是简单线性系数。代码包里一般会有一个read_fy4a.py或类似命名的模块,把h5py.File的读取、通道循环、定标查表封装成函数。你拿到包后第一件事不是跑主程序,而是先确认它用的是h5py还是netCDF4,因为 FY4A 官方分发格式以 HDF5 为主,但有些转存版本会改成 NetCDF,读法完全不同。

2.2 代码包目录的典型构成与职责划分

解压后你大概率会看到类似这样的结构:一个data/或sample/放示例文件,一个src/或根目录放.py脚本,可能还有config.ini或params.py存路径和通道配置。常见的脚本分工是:read_hdf.py负责打开文件和提取数据集,calibration.py负责定标查表,projection.py负责经纬度计算,plot_demo.py负责出图。我拆过不少类似包,最怕的是所有逻辑塞在一个 800 行的脚本里,改一个通道号要翻半天。这份包如果结构清晰,你应该能在calibration.py里找到类似cal_channel(data, cal_table, channel_id)的函数,在projection.py里找到基于NOMSatelliteZenith和NOMSatelliteAzimuth的几何定位代码。先花十分钟把每个文件的头部注释和函数名扫一遍,比直接跑主程序省时间。如果包里没有 README,就看哪个脚本有if __name__ == '__main__',那就是入口。

2.3 依赖环境与版本约束

FY4A 读取代码绕不开的库是h5py、numpy、matplotlib,做投影可能还要pyproj或cartopy。这里有个血泪经验:h5py和numpy的版本兼容性在 Windows 上特别容易翻车,尤其是h5py3.x 配老版本numpy会报DLL load failed。我一般建议用 conda 建独立环境,Python 3.8 到 3.10 之间比较稳。命令如下:

conda create -n fy4a python=3.9 conda activate fy4a conda install h5py numpy matplotlib pyproj cartopy -c conda-forge

逻辑说明:conda-forge渠道的h5py和cartopy预编译包在 Windows 和 Linux 上都比较省心,避免自己编译 GEOS 和 PROJ。参数上,Python 版本不要低于 3.7,因为cartopy新版本已经不支持;也不要盲目上 3.12,部分遥感库还没跟上。装完后用python -c "import h5py; print(h5py.__version__)"验证,能打印版本号再往下走。

3. 读取与定标实战:从 HDF5 到物理量

3.1 用 h5py 提取通道数据与定标表

先看一段典型的读取代码,这段逻辑在多数 FY4A 代码包里都能见到:

import h5py import numpy as np def read_channel(file_path, channel_id): """ 读取指定通道的原始DN值和定标表 channel_id: 1-14, 对应 NOMChannel01-NOMChannel14 """ ch_str = f"{channel_id:02d}" with h5py.File(file_path, 'r') as f: # 原始计数值,通常是uint16 dn = f[f'NOMChannel{ch_str}'][:] # 定标查找表,可见光存反射率,红外存亮温 cal = f[f'CALChannel{ch_str}'][:] # 读取全局属性里的投影参数 attrs = dict(f.attrs) return dn, cal, attrs

逻辑说明:h5py.File以只读模式打开,[:]把数据集完整读进内存。FY4A 全圆盘数据单通道约 2748×2748 像元,uint16 占约 15MB,14 个通道全读约 210MB,内存吃得消。参数上,channel_id用两位字符串格式化是因为 HDF5 里数据集名是NOMChannel01而不是NOMChannel1,少个零就报KeyError。attrs里通常有NOMCenterLat、NOMCenterLon、NOMSatelliteHeight等,后面投影要用。注意有些包会把dn直接叫data,把cal叫lut,看变量名别被绕晕。

3.2 可见光反射率与红外亮温的定标差异

定标是 FY4A 读取里最容易出错的一步。可见光通道的定标表存的是反射率(0-1 或 0-100),红外通道存的是亮温(开尔文)。常见做法是用查找表索引:

def calibrate(dn, cal_table, channel_id): """ 根据通道类型做定标 可见光(1-6): 反射率; 红外(7-14): 亮温 """ dn = dn.astype(np.int32) # 查找表通常是一维数组,长度65536 if channel_id <= 6: result = cal_table[dn] * 100.0 # 转百分比反射率 else: result = cal_table[dn] # 亮温,单位K # 无效值处理,FY4A常用65535或0填充 result[dn >= 65535] = np.nan return result

逻辑说明:dn先转int32是为了避免索引时溢出,uint16最大 65535 刚好是查找表长度上限。可见光乘 100 是习惯做法,有些包不乘,看你后续算法要 0-1 还是 0-100。红外直接查表得亮温,单位是开尔文,要转摄氏度自己减 273.15。无效值判断很关键,FY4A 全圆盘边缘和太空背景常用 65535 填充,不处理的话出图会一片白或一片黑。参数上,cal_table的长度必须是 65536,如果读出来长度不对,说明文件损坏或通道号错了。

3.3 批量处理多个文件的循环写法

单文件跑通后,实际业务都是批量。常见做法是用glob或os.listdir遍历:

import glob import os def batch_process(data_dir, channel_id, out_dir): files = sorted(glob.glob(os.path.join(data_dir, '*.HDF'))) for fp in files: dn, cal, attrs = read_channel(fp, channel_id) result = calibrate(dn, cal, channel_id) # 用文件名做输出名,避免覆盖 base = os.path.splitext(os.path.basename(fp))[0] np.save(os.path.join(out_dir, f'{base}_ch{channel_id:02d}.npy'), result)

逻辑说明:glob匹配.HDF后缀,FY4A 官方文件后缀有时是大写.HDF,有时是.hdf,Windows 不区分但 Linux 区分,建议用*.HDF*或统一转小写。sorted保证时间顺序,做时间序列时很重要。存.npy比存图片快,后续要用再读。参数上,out_dir要先os.makedirs(out_dir, exist_ok=True),不然第一次跑就报目录不存在。如果文件多,可以加tqdm看进度,但别在正式流程里留 print,日志用logging。

4. 经纬度定位与投影:把像元对上地图

4.1 基于卫星天顶角和方位角的几何定位

FY4A 的定位信息藏在NOMSatelliteZenith和NOMSatelliteAzimuth两个数据集里,有些包还会用NOMCenterLat、NOMCenterLon和NOMSatelliteHeight做球面三角计算。常见做法是直接读角度数据集,再用简单几何关系反算经纬度:

def compute_latlon(zenith, azimuth, sat_height=35786.0, center_lat=0.0, center_lon=104.7): """ 简化几何定位,适用于全圆盘 zenith/azimuth: 角度,单位度 sat_height: 卫星高度,单位km """ import numpy as np re = 6371.0 # 地球半径km zen_rad = np.deg2rad(zenith) azi_rad = np.deg2rad(azimuth) # 星下点距离角 gamma = np.arcsin((re + sat_height) / re * np.sin(zen_rad)) - zen_rad lat = np.arcsin(np.sin(center_lat) * np.cos(gamma) + np.cos(center_lat) * np.sin(gamma) * np.cos(azi_rad)) lon = center_lon + np.arctan2(np.sin(azi_rad) * np.sin(gamma) * np.cos(center_lat), np.cos(gamma) - np.sin(center_lat) * np.sin(lat)) return np.rad2deg(lat), np.rad2deg(lon)

逻辑说明:这是简化球面模型,精度对全圆盘预览够用,做定量分析建议用官方投影库。sat_height取 35786km 是 FY4A 标称高度,center_lon取 104.7 是 FY4A 星下点经度。参数上,zenith和azimuth数据集可能叫NOMSatelliteZenith或NOMSatelliteZenithAngle,看文件里实际 key。如果包里用的是pyproj的Geostationary投影,那更准,但配置参数多,新手容易配错。

4.2 用 cartopy 出图的参数配置

出图是验证读取是否正确的最快方式。用cartopy画全圆盘:

import matplotlib.pyplot as plt import cartopy.crs as ccrs def plot_disk(data, lat, lon, out_png): fig = plt.figure(figsize=(10, 10)) ax = plt.axes(projection=ccrs.PlateCarree()) ax.coastlines(resolution='110m') # 注意经纬度网格要和数据形状一致 im = ax.pcolormesh(lon, lat, data, transform=ccrs.PlateCarree(), cmap='gray', vmin=0, vmax=100) plt.colorbar(im, ax=ax, shrink=0.7) plt.savefig(out_png, dpi=150, bbox_inches='tight') plt.close()

逻辑说明:PlateCarree是等经纬度投影,适合快速预览,但全圆盘边缘会变形。pcolormesh比imshow更适合不规则网格,但要求lon、lat、data形状一致。参数上,vmin、vmax按通道调,可见光 0-100,红外 180-320。coastlines的resolution用110m最快,50m更细但慢。如果出图空白,先检查data是不是全nan,再检查lon、lat有没有正确从角度算出来。

4.3 投影参数表与常见取值

参数典型值说明
卫星高度35786 kmFY4A 标称
星下点经度104.7°E全圆盘中心
地球半径6371 km简化球模型
通道空间分辨率0.5-4 km可见光高,红外低
无效值65535uint16 填充

这张表建议存成config.ini,别硬编码在脚本里。换星或换区域产品时只改配置,不动代码。

5. 避坑与排查:那些让代码跑不通的细节

5.1 现象:打开文件报 KeyError: 'NOMChannel01'

原因:FY4A 文件里数据集名可能带前缀或后缀,比如NOMChannel01写成NOMChannel_01,或者区域产品用Channel01。解决:先用h5py.File(fp,'r').keys()打印所有 key,确认实际命名,再改代码里的字符串拼接。别猜,直接看。

5.2 现象:定标后数据全是 0 或全是 65535

原因:查找表索引越界或无效值判断写反。cal_table[dn]如果dn里有大于表长度的值会报错,但有些包用np.clip把越界值截断,结果全落在无效区。解决:先检查dn.max()和cal_table.shape,确认dn最大值小于表长度。无效值判断用dn == 65535而不是dn >= 65535,因为有些文件用 0 填充。

5.3 现象:出图只有一小块有数据,其余全黑

原因:经纬度计算时用了错误的星下点经度或卫星高度,导致投影偏移。解决:确认center_lon是 104.7 而不是 0,sat_height是 35786 而不是 42164(那是别的星)。如果包里用pyproj,检查Geostationary的satellite_height参数单位是米还是千米。

5.4 现象:批量处理到一半内存溢出

原因:循环里没释放大数组,或者h5py文件句柄没关。解决:用with h5py.File(...) as f确保关闭,循环里del dn, cal, result手动释放,或者用gc.collect()。全圆盘 14 通道全读再处理,内存峰值可能上 GB,分批处理更稳。

5.5 现象:cartopy 出图报 PROJ_LIB 错误

原因:conda 环境里pyproj和cartopy的 PROJ 数据路径冲突。解决:在脚本开头加import os; os.environ['PROJ_LIB'] = '/path/to/conda/env/share/proj',路径用conda env list找到环境目录后拼。或者干脆用pyproj自己算经纬度,出图用matplotlib原生imshow,绕开cartopy。

6. 进阶技巧:把读取代码接进自己的流程

跑通单文件后,下一步是把它变成可复用的模块。我一般会把read_channel、calibrate、compute_latlon三个函数抽到一个fy4a_reader.py里,用argparse做命令行入口,这样既能被其他脚本import,也能直接跑。参数配置全部走config.ini,通道号、输入目录、输出目录、投影参数都放进去,换数据不用改代码。验证方法很简单:拿一个已知时间的全圆盘文件,读通道 12(红外窗区),定标后看亮温范围是不是在 180-320K 之间,出图看青藏高原是不是比同纬度海洋冷。如果对不上,先查定标表,再查经纬度。还有一个习惯:每次处理新一批数据前,先拿一个文件跑h5py的keys()和attrs打印,确认格式没变。FY4A 数据版本迭代时,数据集命名和属性偶尔会调整,不先检查就批量跑,跑完才发现全错,后悔药没处买。从那以后我每次接新数据源,都强制走一遍“打印 keys → 读单通道 → 定标 → 出图 → 对已知地物”这五步,确认无误再上批量。希望帮到你。

本文还有配套的精品资源,点击获取

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询