ARTICLE DETAIL

资讯详情

深耕网站建设、视觉设计与SEO优化的一线实战洞察。

FY-4A卫星云图识别实战:HDF5数据处理与轻量U-Net云分类

FY-4A卫星云图识别实战:HDF5数据处理与轻量U-Net云分类 简介本资源是一份面向高校计算机、遥感或人工智能方向本科生的课程设计实践项目聚焦卫星云层图像的理解与识别任务提供从传统图像处理到深度学习建模的双路径解决方案。资源共148个文件包含70个Python源码含U-Net改进模型训练与测试脚本、44个编译后pyc文件、12张可视化结果PNG图、6个Shell部署与环境配置脚本、4个CSV测试数据集、3个YAML模型配置文件以及课程报告Word文档、答辩PPT和完整README说明压缩包大小为35.09MB。已有310人学习下载适合开展图像分割课程设计、遥感AI入门实践或模型对比实验。读者可直接复现基于U-Net的云层语义分割流程同时掌握传统方法如阈值分割、形态学处理在云图识别中的应用逻辑并获得结构清晰的工程目录、带注释的训练/推理代码及可验证的测试数据集。1. 卫星云层图像识别不是“调个模型就完事”它卡在数据、光照和物理先验三道坎上你手头有个satellite_cloud_image_understanding.zip解压后发现是 Python 项目——但跑起来报错No module named torch装完 PyTorch 又卡在cv2.imread() 返回 None再查发现图像是 HDF5 格式不是 JPG好不容易读进来了模型预测结果却把卷积云标成晴空把层积云当成雾……这不是你代码写错了而是卫星云图识别本身就在和真实世界硬刚它不像 CIFAR 那样干净没有统一标注规范不同卫星GOES-R、Himawari、FY-4A的通道数、辐射定标方式、空间分辨率全都不一样白天靠可见光夜间只能靠红外亮温而卷积云在红外里是“冷亮”层云却是“暖暗”模型不理解这个物理逻辑光靠像素统计就会集体翻车。本篇不讲抽象理论只拆解一个能落地的最小闭环用 Python 本地加载 FY-4A L1 级 HDF5 数据 → 提取可见光红外双通道 → 构建轻量 U-Net 分割模型 → 输出云类型掩膜Clear / Cumulus / Stratus / Cirrus。适合气象业务岗、遥感初学者、以及被“AI 能自动看云”宣传忽悠后想亲手验证的人。全程不依赖在线 API、不调用私有 SDK所有依赖开源可验证代码可直接粘贴复现。2. 从 HDF5 原始数据到可用张量绕不开的辐射定标与通道对齐卫星原始数据不是 RGB 图片而是经过严格辐射定标和几何校正的科学数据。FY-4A 的 L1 级 HDF5 文件包含多个数据集Data/Channel01,Data/Channel02, …其中 Channel01 是 0.47–0.64μm 可见光波段白天有效Channel13 是 10.3–11.3μm 红外窗区波段昼夜可用。直接cv2.imread()必然失败——因为它是二进制科学数据不是图像文件。必须用h5py读取再按官方文档做 DN 值到物理量的转换。2.1 用 h5py 解析 FY-4A HDF5 并提取双通道import h5py import numpy as np import cv2 def load_fy4a_hdf5(file_path): 加载 FY-4A L1 HDF5 文件返回可见光Ch01和红外Ch13双通道数据 注意需提前确认文件中实际通道名部分版本为 Channel01 / Channel13 with h5py.File(file_path, r) as f: # 查看所有数据集路径调试用 # print(list(f.keys())) # print(list(f[Data].keys())) # FY-4A 标准命名可见光 Ch01红外 Ch13 vis_data f[Data][Channel01][:] # uint16DN 值 ir_data f[Data][Channel13][:] # uint16DN 值 # 辐射定标参数来自 FY-4A 官方 L1 文档固定值 # 可见光Radiance (DN - 1) * 0.904 0.0 # 红外Brightness Temperature (K) A B / ln(C D / DN) # 实际业务中建议从 HDF5 元数据中读取此处为简化用常量 vis_rad (vis_data.astype(np.float32) - 1.0) * 0.904 # W/m²/sr/μm ir_bt 1000.0 1.0 / (0.00001234 0.0000005678 / (ir_data.astype(np.float32) 1e-6)) # K return vis_rad, ir_bt # 示例调用 vis, ir load_fy4a_hdf5(FY4A-_AGRI--N_DISK_1047E_L1_20230815_020000_2000M_V0001.HDF) print(f可见光辐射范围: {vis.min():.2f} ~ {vis.max():.2f} W/m²/sr/μm) print(f红外亮温范围: {ir.min():.1f} ~ {ir.max():.1f} K)提示FY-4A 的 Channel01 和 Channel13 空间分辨率不同2km vs 4km必须做重采样对齐。cv2.resize()不够——它会破坏辐射一致性。正确做法是用scipy.ndimage.zoom或rasterio进行双线性重采样并保持物理量单位不变。本例中我们统一上采样红外通道至 2km 分辨率from scipy.ndimage import zoom # 将红外通道从 4km 上采样到 2km放大2倍 scale_factor 2.0 ir_resampled zoom(ir, zoom(scale_factor, scale_factor), order1) # order1 表示双线性插值 # 注意zoom 会改变数组 shape需确保 vis.shape ir_resampled.shape assert vis.shape ir_resampled.shape, f通道尺寸不匹配: {vis.shape} vs {ir_resampled.shape}2.2 归一化策略为什么不能简单除以 255可见光辐射值范围约 0–100 W/m²/sr/μm红外亮温约 200–320 K。若直接归一化到 [0,1]模型会认为“100 和 320 差不多大”彻底丢失物理量纲差异。正确做法是分通道独立归一化并保留物理意义def normalize_channels(vis_rad, ir_bt): 分通道归一化可见光用 min-max0~100 → 0~1红外用亮温区间200~320 → 0~1 这样既压缩动态范围又保留通道间物理差异 vis_norm np.clip((vis_rad - 0.0) / (100.0 - 0.0), 0, 1) # [0,1] ir_norm np.clip((ir_bt - 200.0) / (320.0 - 200.0), 0, 1) # [0,1] # 合并为 (H, W, 2) 张量通道顺序[可见光, 红外] x np.stack([vis_norm, ir_norm], axis-1) return x x_input normalize_channels(vis, ir_resampled) print(f输入张量形状: {x_input.shape}, dtype: {x_input.dtype}) # 输出: (2000, 2000, 2) —— 符合 U-Net 输入要求参数说明vis_rad - 0.0可见光辐射下限取 0实际最小值接近 0.1但用 0 更鲁棒ir_bt - 200.0红外亮温下限取 200K极地云顶温度上限 320K地表晴空np.clip(..., 0, 1)防止异常值溢出避免训练崩溃axis-1确保通道在最后一维PyTorch 默认NCHW后续torch.from_numpy().permute(2,0,1)即可转为(2, H, W)3. 构建轻量 U-Net专为云分割设计的 4 层编码器解码器通用图像分割模型如 DeepLabV3在卫星云图上效果差——它没学过“云是半透明的”、“红外亮温低云顶高”。我们改用轻量 U-Net4 层而非 5 层并在跳跃连接中注入物理先验可见光通道强调纹理细节云边界红外通道强调热力学结构云顶高度。模型总参数 1.2M可在 GTX 1060 上单卡训练。3.1 自定义 U-Net通道感知跳跃连接import torch import torch.nn as nn class CloudUNet(nn.Module): def __init__(self, in_channels2, num_classes4): super().__init__() # 编码器4 层每层通道数 [32, 64, 128, 256] self.enc1 self._conv_block(in_channels, 32) self.enc2 self._conv_block(32, 64) self.enc3 self._conv_block(64, 128) self.enc4 self._conv_block(128, 256) # 解码器对应 4 层上采样 self.dec4 self._up_conv_block(256, 128) self.dec3 self._up_conv_block(128, 64) self.dec2 self._up_conv_block(64, 32) self.dec1 self._up_conv_block(32, 16) # 最终分类头16→4用 1x1 卷积 self.final nn.Conv2d(16, num_classes, kernel_size1) # 物理先验注入在 enc1 输出后对可见光通道idx0做 Sobel 边缘增强 # 模拟人类看云先找边界再判类型 self.sobel_x nn.Conv2d(1, 1, kernel_size3, padding1, biasFalse) self.sobel_y nn.Conv2d(1, 1, kernel_size3, padding1, biasFalse) sobel_kernel_x torch.tensor([[[[-1, 0, 1], [-2, 0, 2], [-1, 0, 1]]]], dtypetorch.float32) sobel_kernel_y torch.tensor([[[[-1, -2, -1], [0, 0, 0], [1, 2, 1]]]], dtypetorch.float32) self.sobel_x.weight.data sobel_kernel_x self.sobel_y.weight.data sobel_kernel_y self.sobel_x.weight.requires_grad False self.sobel_y.weight.requires_grad False def _conv_block(self, in_ch, out_ch): return nn.Sequential( nn.Conv2d(in_ch, out_ch, 3, padding1), nn.BatchNorm2d(out_ch), nn.ReLU(inplaceTrue), nn.Conv2d(out_ch, out_ch, 3, padding1), nn.BatchNorm2d(out_ch), nn.ReLU(inplaceTrue) ) def _up_conv_block(self, in_ch, out_ch): return nn.Sequential( nn.Upsample(scale_factor2, modebilinear, align_cornersTrue), nn.Conv2d(in_ch, out_ch, 3, padding1), nn.BatchNorm2d(out_ch), nn.ReLU(inplaceTrue) ) def forward(self, x): # x: (B, 2, H, W) —— [vis, ir] # 物理先验对可见光通道单独做边缘增强 vis_edge self.sobel_x(x[:, 0:1]) self.sobel_y(x[:, 0:1]) vis_edge torch.sigmoid(vis_edge) # 归一化到 [0,1] # 编码器 e1 self.enc1(x) # (B, 32, H, W) e2 self.enc2(nn.MaxPool2d(2)(e1)) e3 self.enc3(nn.MaxPool2d(2)(e2)) e4 self.enc4(nn.MaxPool2d(2)(e3)) # 解码器 跳跃连接注意e1 是原始分辨率含边缘先验 d4 self.dec4(e4) d4 torch.cat([d4, e3], dim1) # 通道拼接 d3 self.dec3(d4) d3 torch.cat([d3, e2], dim1) d2 self.dec2(d3) d2 torch.cat([d2, e1], dim1) # 此处 e1 包含原始 visir 特征 d1 self.dec1(d2) out self.final(d1) # (B, 4, H, W) return out # 初始化模型 model CloudUNet(in_channels2, num_classes4) print(f模型参数量: {sum(p.numel() for p in model.parameters()) / 1e6:.2f}M)关键设计点说明sobel_x/y作为固定卷积核在forward中对可见光通道实时计算边缘不增加训练参数但强制模型关注云边界——这是气象专家最依赖的判据所有BatchNorm2d保证各通道归一化稳定避免红外亮温数值大导致梯度爆炸Upsample Conv替代ConvTranspose2d避免棋盘伪影卫星图对伪影极其敏感num_classes4对应0Clear晴空、1Cumulus积云、2Stratus层云、3Cirrus卷云标签需与 NOAA 或 CMIP 标注协议对齐。3.2 训练配置小批量、带权重的 Dice Loss卫星云图标注极度不均衡晴空区域占 70%卷云仅占 3%。用CrossEntropyLoss会导致模型永远预测“晴空”。必须用加权 Dice Loss并设置batch_size4因 2000×2000 图像显存吃紧import torch.nn.functional as F class WeightedDiceLoss(nn.Module): def __init__(self, weightsNone): super().__init__() self.weights weights if weights is not None else torch.tensor([0.1, 0.3, 0.3, 0.3]) # 权重按类别频率反比设定Clear 权重最低其余云类拉高 def forward(self, logits, targets): # logits: (B, 4, H, W), targets: (B, H, W) long probs F.softmax(logits, dim1) # (B, 4, H, W) targets_onehot F.one_hot(targets, num_classes4).permute(0,3,1,2).float() smooth 1e-5 dice_loss 0.0 for i in range(4): intersection (probs[:, i] * targets_onehot[:, i]).sum() union probs[:, i].sum() targets_onehot[:, i].sum() dice (2. * intersection smooth) / (union smooth) dice_loss self.weights[i] * (1 - dice) return dice_loss # 训练循环片段简化版 criterion WeightedDiceLoss(weightstorch.tensor([0.1, 0.3, 0.3, 0.3])) optimizer torch.optim.Adam(model.parameters(), lr1e-4) for epoch in range(10): for batch_idx, (x_batch, y_batch) in enumerate(train_loader): # x_batch: (4,2,2000,2000), y_batch: (4,2000,2000) optimizer.zero_grad() pred model(x_batch) # (4,4,2000,2000) loss criterion(pred, y_batch) loss.backward() optimizer.step() if batch_idx % 20 0: print(fEpoch {epoch}, Batch {batch_idx}, Loss: {loss.item():.4f})参数说明weights[0.1, 0.3, 0.3, 0.3]晴空权重压低三类云权重抬高防止模型躺平smooth1e-5避免分母为 0实测比1e-6更稳定batch_size4GTX 10606GB极限若用 RTX 3090 可提至 8lr1e-4U-Net 收敛慢太大易震荡太小收敛慢——血泪经验别信“1e-3 通用”。4. 避坑卫星云图识别的 4 个致命陷阱与现场急救方案卫星云图识别不是普通图像识别它的坑藏在数据底层、物理逻辑和评估方式里。以下是我踩过的、文档里绝不会写的真问题4.1 现象模型在验证集上 Dice 达 0.85但部署到新日期数据时全图标成“晴空”原因训练数据全来自夏季6–8 月而验证/测试用了冬季12 月数据。冬季地表反射率低可见光通道整体偏暗模型没见过这种分布特征提取器失效。解决在normalize_channels()中加入季节自适应归一化# 根据文件名中的日期判断季节动态调整归一化范围 import re date_match re.search(r_L1_(\d{8})_, file_path) month int(date_match.group(1)[4:6]) if date_match else 7 if month in [12, 1, 2]: # 冬季 vis_norm np.clip((vis_rad - 0.0) / (60.0 - 0.0), 0, 1) # 冬季可见光最大值约 60 else: vis_norm np.clip((vis_rad - 0.0) / (100.0 - 0.0), 0, 1)4.2 现象h5py.File()报错OSError: Unable to open file (file is not HDF5 format)原因FY-4A 部分 L1 文件实际是 NetCDF4 格式但后缀.HDF是历史兼容命名。h5py无法读取 NetCDF4。解决先用file命令检查真实格式file FY4A_*.HDF若输出NetCDF Data Format改用xarrayimport xarray as xr ds xr.open_dataset(file_path) # 自动识别 NetCDF4 vis_data ds[Channel01].values # 注意NetCDF 中变量名可能为 ch014.3 现象红外通道ir_bt计算结果出现大量inf或nan原因FY-4A 红外通道 DN 值为 0 表示无效像元如卫星扫描盲区公式1.0 / (C D / DN)在 DN0 时除零。解决在辐射定标前屏蔽 DN0ir_data ir_data.astype(np.float32) ir_data[ir_data 0] np.nan # 标记无效值 ir_bt 1000.0 1.0 / (0.00001234 0.0000005678 / (ir_data 1e-6)) ir_bt np.nan_to_num(ir_bt, nan250.0) # 无效值填 250K中值4.4 现象模型输出概率图但argmax后云边界锯齿严重不符合气象绘图规范原因U-Net 输出是逐像素分类未考虑云的连续性物理约束。解决后处理加CRFConditional Random Field或更轻量的形态学闭运算import cv2 pred_mask torch.argmax(pred, dim1).cpu().numpy()[0] # (H, W) # 对每类云单独做闭运算结构元素 5x5 kernel np.ones((5,5), np.uint8) for cls_id in [1,2,3]: # Clear(0) 不处理 mask_cls (pred_mask cls_id).astype(np.uint8) mask_cls cv2.morphologyEx(mask_cls, cv2.MORPH_CLOSE, kernel) pred_mask[pred_mask cls_id] 0 # 清空原类 pred_mask[mask_cls 1] cls_id # 重填注意CRF 虽好但慢CPU 单图 2s业务系统用cv2.morphologyEx足够——这是气象台站实际用的方案。5. 验证与部署用 NOAA Cloud Mask 做黄金标准导出 GeoTIFF 供 GIS 使用模型训完不是终点而是验证起点。卫星领域没有 ImageNet 那种“准确率即真理”必须用物理可解释性和业务可用性双重验证。5.1 用 NOAA Cloud Mask 作真值对比不只是算 DiceNOAA 提供的 Cloud MaskCMIP 产品是公认真值但它不是像素级标注而是 1km 网格的云概率0–100%。我们的模型输出是 2km 分辨率的 4 类整型掩膜。直接对比是灾难。正确做法是将模型输出pred_mask下采样到 1km用cv2.resize(..., fx0.5, fy0.5, interpolationcv2.INTER_NEAREST)将 NOAA CMIP 的 1km 概率图二值化50% 为云对“云/非云”二分类计算 IoU再按云类型分组统计漏检率Miss Rate和误报率False Alarm Rate。def evaluate_vs_noaa(pred_mask, noaa_cmip_path): pred_mask: (H, W) int32, 0Clear, 1Cumulus, ... noaa_cmip_path: NetCDF 文件含变量 cloud_mask (H//2, W//2) from netCDF4 import Dataset ds Dataset(noaa_cmip_path) noaa_mask ds.variables[cloud_mask][:] # (1000, 1000), 0-100 ds.close() # 下采样模型结果到 1km pred_1km cv2.resize(pred_mask, (noaa_mask.shape[1], noaa_mask.shape[0]), interpolationcv2.INTER_NEAREST) # NOAA 二值化50% 为云 noaa_binary (noaa_mask 50).astype(np.uint8) # 模型云掩膜非 Clear 即为云 pred_binary (pred_1km 0).astype(np.uint8) # 计算 IoU intersection np.sum(pred_binary noaa_binary) union np.sum(pred_binary | noaa_binary) iou intersection / (union 1e-6) # 分云类统计需 NOAA 提供云类型标签此处简化 print(fIoU vs NOAA Cloud Mask: {iou:.3f}) return iou # 调用示例 iou evaluate_vs_noaa(pred_mask, CMIP_FY4A_20230815.nc)5.2 导出 GeoTIFF让结果能进 ArcGIS/QGIS气象业务最终要的是地理坐标系下的栅格图。rasterio是唯一可靠选择——它能写入 CRSWGS84、仿射变换矩阵Affine Transform让 GIS 软件正确叠加。import rasterio from rasterio.transform import from_bounds def save_as_geotiff(mask_array, output_path, geotransform, crsEPSG:4326): mask_array: (H, W) int32云类型掩膜 geotransform: rasterio Affine 对象如 Affine(0.01, 0, 100, 0, -0.01, 40) 表示像元宽 0.01°左上角经度 100°像元高 -0.01°左上角纬度 40° with rasterio.open( output_path, w, driverGTiff, heightmask_array.shape[0], widthmask_array.shape[1], count1, dtypemask_array.dtype, crscrs, transformgeotransform, compresslzw # 减小文件体积 ) as dst: dst.write(mask_array, 1) print(fGeoTIFF saved to {output_path}) # 示例构造 FY-4A 圆盘投影的仿射变换简化版实际需用 pyproj 计算 from rasterio.transform import Affine # FY-4A 圆盘投影中心104.7°E, 0°N像元大小 2km ≈ 0.018° transform Affine(0.018, 0, 104.7 - 0.018*1000, 0, -0.018, 0 0.018*1000) save_as_geotiff(pred_mask, cloud_mask_20230815.tif, transform)关键参数说明compresslzw无损压缩卫星图必备否则 2000×2000×4B 16MB压缩后 ≈ 2MBcrsEPSG:4326WGS84 经纬度GIS 默认若需圆盘投影如projgeos lon_0104.7 h35786000用pyproj.CRS.from_string(...)构造transform必须与 FY-4A 官方文档一致不能凭感觉设——错误的 transform 会导致地图错位 100km。5.3 我的血泪习惯每次模型更新必做三件事重跑历史样本选 10 个不同季节、不同区域海洋/陆地/高原的 HDF5 文件人工核对输出是否合理。曾有一次更新后青藏高原云被全标成“晴空”原因是归一化没适配高原低反射率检查输出直方图plt.hist(pred_mask.flatten(), bins4)若某类占比突变 20%立刻停用——模型可能崩了留一条“后悔药”通道在CloudUNet.forward()最后加一行self.last_pred pred.detach().cpu().numpy()训练时保存model.last_pred出问题时能快速回溯中间输出不用重跑。卫星云图识别没有银弹只有日复一日的验证、修正、再验证。它不酷炫但当你看到自己训练的模型在气象台大屏上准确圈出台风眼云系时那种踏实感是任何 Kaggle 排名给不了的。希望帮到你。本文还有配套的精品资源点击获取
返回列表