企业官网建设流程全解析

简介：面向气象领域与雷达数据处理学习者的多普勒雷达强度显示工具，基于C++工程实现，针对敏视达气象雷达的回波数据，覆盖从脉冲信号采集、距离解析、多普勒速度计算、反射率强度分析到彩色雷达图生成的完整流程。压缩包共13个文件、大小约737KB，构成一个完整的C++雷达数据处理工程；内容涵盖dat雷达基数据、cpp源代码、exe可执行程序、dsw/dsp工程配置以及obj、pdb等构建与调试辅助文件，既能直接运行查看雷达强度分布，也便于在Visual C++环境中对照源码二次开发。目前已有258人学习下载。包内附带的多个实测数据文件，可帮助读者结合源码理解降水粒子速度反演、强度分级定标与雷达图着色原理，适合气象专业学生、科研人员及雷达数据爱好者用于教学演示、算法验证和科研分析。

1. 多普勒气象雷达强度数据，为什么ShowRadarData_intensity值得自己写一套

拿到多普勒气象雷达基数据之后，第一件让人头疼的事往往不是算法，而是数据根本打不开。厂家的商业软件把雷达的每一个仰角、每一条径向数据都封装在黑匣子里，你只知道屏幕上有一张漂亮的回波图，却不知道强度值是怎么从文件里来的。ShowRadarData_intensity这类以强度场为核心的脚本，就是为了解决这个“从原始数据到气象雷达平面显示ppi”的断层。我第一次用类似脚本跑出第一张强度图时，整片螺旋回波从雷达中心向外展开，那一刻才真正敢说自己在处理气象雷达数据。这篇笔记适合刚拿到雷达基数据、想绕开商业软件自己画PPI的数据工程师，也适合需要批量出图做质控的预报支撑岗位。

2. 读雷达基数据之前：先搞懂格式和强度字段

2.1 三种常见基数据格式：SXD、CINRAD、NetCDF，认错等于白干

多普勒雷达基数据的格式从来不是一个通用标准，我经手过的文件至少能分成三类：

一是老一代SXD格式，常见于部分早期多普勒雷达，整份体扫数据按仰角顺序连续存放，每个仰角块里有若干条径向，径向之间用固定长度的头部字段分隔。这种格式没有明显的文件头字符串，需要根据帧同步码找起始位置，偏移量算错一位后面全是乱码。

二是新一代CINRAD格式，国内SA/SB型业务雷达用的就是这一类，文件以二进制分块存储，径向头包含仰角、方位角、库数、脉冲数等信息，身体里再按强度、速度、谱宽三个数据块排列。这个格式相对规范，但版本多，新旧雷达的径向头长度不一样，读取时最怕按旧版本写死长度。

三是NetCDF格式，现在不少新型相控阵和天气雷达直接输出标准NetCDF，变量名和维度一眼就能看明白，比如reflectivity、azimuth、elevation、range。处理这种格式最省心，但问题在于字段名并不统一，有的写REF，有的写Z_HDR。

拿到一个陌生文件，我一般先用ncdump -h看一遍，如果是NetCDF直接列出了维度和变量就很幸运；如果是二进制文件，就用xxd看前256字节，找ASCII字符串和固定头长度。这一步必须花十分钟把格式认准，否则后面画出来的图可能就是噪点，到时候还得回过来查格式，浪费时间。ShowRadarData_intensity这类脚本里通常内置了若干解析器，但你要知道它在解哪一种格式，才能把参数配到对应的分支上。

2.2 intensity到底是哪个字段：反射率因子和强度码的换算

气象雷达平面显示ppi上的颜色，核心就是“强度”这个物理量，但文件里的intensity字段往往不是直接可用的dBZ。我见过太多次别人把强度码当成强度值直接画图，出来的图颜色完全错乱。

常见存储方式是：每一个距离库上存一个字节或两个字节的强度码，取值从0到255，真正气象反射率因子dBZ需要做线性换算。典型换算公式是dBZ = value * 0.5 - 33.0，也就是说值0对应-33dBZ，值66对应0dBZ，值100对应17dBZ。但这个斜率0.5和截距-33不是雷打不动的，不同厂家的雷达在标定时会调截距，我见过-32、-34，甚至用0.01和-0.01的微调版本。ShowRadarData_intensity如果提供配置文件，通常会有slope和intercept这两个参数，必须根据自己的雷达标定值改。

还有更隐蔽的坑：强度码里有保留值。很多格式规定1代表缺测，2代表地物杂波，0是无效数据，这些不能参与线性换算，否则画出来在零度层附近会出现一堆莫须有的-33dBZ像素。真正有效的降水回波强度码大多在66到180之间，换算后就是0到57dBZ。如果你统计出来的强度码分布大量集中在200以上或30以下，先别画图，回去检查字节序和符号位，大概率读错了位置。

2.3 解析骨架：按字节拆径向数据，先跑通一个仰角

写解析代码不需要把整个格式规范背下来，只需要拆出一个仰角的数据就能验证思路。下面是个简化版CINRAD格式读取骨架，演示怎么从径向头走到强度数据区。

import struct import numpy as np def read_radial(fp): # 径向头：固定128字节，版本不同可能变化 head_bytes = fp.read(128) if len(head_bytes) < 128: return None # 按小端序解析前12字节 elevation, azimuth, n_gate = struct.unpack('<ffI', head_bytes[:12]) # 读取n_gate个强度码，按uint16读取 code = np.frombuffer(fp.read(n_gate * 2), dtype='<u2').astype(float) # 做缺测掩膜：值为1或2的剔除 code[(code == 1) | (code == 2)] = np.nan # 强度码转dBZ，斜率和截距从配置文件读取 dbz = code * 0.5 - 33.0 return elevation, azimuth, dbz

逻辑说明：先把径向头里仰角、方位角、库数拆出来，再按库数读取相同长度的强度数据。代码里用np.nan覆盖特殊标记，避免绘图时把杂波当降水画出来。字节序要跟文件一致，'<'是小端，如果文件是大端就得换成'>'，这一项不对直接把所有值读成乱数据。

参数说明：n_gate是距离库数，通常一条径向有300到1000个库；elevation和azimuth单位是度，后面投影前要转弧度。斜率和截距建议从外部参数文件读，不要写在函数里。这段逻辑跑同一个仰角的任意一条径向，如果算出来的最大值在50到60左右，且最小值不是极端负数，基本就对了。

3. 用ShowRadarData_intensity画出第一张PPI：从解压到出图

3.1 解压之后先别急着跑：看脚本入口和依赖

网上常见的ShowRadarData_intensity压缩包解压之后，里面通常是一个Python项目，包含read_radar.py、plot.py、config.ini和示例数据目录。我拿到手不会直接python乱跑，先打开README或main函数看依赖。它一般需要numpy、matplotlib、scipy，如果画地图可能还需要cartopy或basemap。缺少依赖时直接用会报ModuleNotFoundError，所以先用pip list确认环境。

另外检查一下脚本有没有硬编码路径。很多脚本在设计时会把输入文件路径写在代码开头，比如file_path = 'D:\\Radar\\20240815.dat'，这种情况下先把绝对路径改成相对路径或者命令行参数。我习惯的做法是给脚本套一层薄薄的命令行入口，不直接改原始脚本，方便以后换目录结构。

# 先查看帮助 python show_intensity.py --help # 常见参数形式 python show_intensity.py --input radar.dat --output ppi.png --radius 150

如果脚本完全没有命令行接口，就自己写一个run.py去调用它的函数，这样保留原始脚本的纯净性，也方便版本升级后对比差异。

3.2 最小可运行命令和函数调用

无论脚本包装成什么样，核心逻辑基本是“读取基数据，提取某一仰角的强度场，然后绘制极坐标投影图”。我用过的最简命令是把读取和绘图分开：

from show_intensity import load_radar_file, plot_ppi # 读取整个体扫，返回数据体和元信息 data, meta = load_radar_file('20240815_0200.dat') # 把第1个仰角（通常是最低仰角）画成PPI plot_ppi(data[0], meta, 'ppi_e001.png')

逻辑说明：load_radar_file内部会解析文件里所有仰角，返回一个三维数组，第一维是仰角索引，第二维是径向方位，第三维是距离库。meta里保存雷达经纬度、海拔、库长、仰角列表等关键元数据。plot_ppi负责把强度数组按极坐标投影到平面，并套用色标。

参数说明：如果只关心某一个仰角，比如0.5度，索引一般就是0，但有的文件最低仰角是0.0度，要看meta['elevations']确认。plot_ppi内部通常会自动取最大距离范围，如果没有参数控制，默认画满整个文件量程，但这样边缘回波会被压缩，建议手动限制半径。

3.3 颜色条和坐标投影：让气象雷达平面显示ppi像业务软件

PPI的“像不像气象产品”，一半靠色标。matplotlib默认的jet色标把低值画成蓝色高值画成红色，但气象上反射率因子有固定约定：负值到0dBZ用蓝绿，10到30dBZ用绿色，35以上用黄橙红，超过50用紫色或洋红。很多业务软件的颜色表是自定义的，ShowRadarData_intensity如果没有附带色标文件，就需要自己做一个离散色标。

import matplotlib.colors as mcolors levels = [-20, 0, 5, 10, 15, 20, 25, 30, 35, 40, 45, 50, 55, 60, 65, 70] colors = [ '#0055CC', '#0088FF', '#00BBBB', '#00DD44', '#22BB22', '#88CC22', '#CCCC00', '#FFAA00', '#FF6600', '#FF3300', '#CC0000', '#9900CC', '#6600AA', '#330077' ] cmap = mcolors.ListedColormap(colors) norm = mcolors.BoundaryNorm(levels, cmap.N)

逻辑说明：这里把色标分成16个级别，每个区间一个颜色，比连续色标更符合业务习惯。使用BoundaryNorm会将数据值映射到对应区间，不足-20dBZ和超过70dBZ的像素统一变成边界色。绘图时把cmap和norm传给pcolormesh，色标颜色就不会跨度太大。

坐标投影方面，极坐标转平面是PPI的关键。通常在笛卡尔坐标系里创建两个网格，然后算每个网格对应的距离和方位角，再从强度数组中插值。反插值方法简单但会损失细节，常见做法是用scipy.interpolate.griddata把极坐标点插到规则网格上。网格分辨率一般设为雷达库长的五分之一到十分之一，比如库长1km，网格设200m到500m，太高反而插值出很多空洞。

3.4 三个必调参数：雷达中心、最大距离、dBZ范围

先看雷达中心。如果脚本读出来的是极坐标相对位置，绘图时得把雷达经纬度和地图投影绑定，否则叠加地图时回波错位能差几十公里。我在配置里固定写radar_lat和radar_lon，同时确认单位是十进制度还是度分秒格式。一个常见低级错误是把经纬度写反，北半球纬度范围90到0，经度范围180到-180，如果纬度小于20且经度大于100，多半写反了。

再看最大距离。雷达体扫通常能探测到230km甚至300km，但PPI显示范围超过200km后，低仰角波束会被地球曲率抬高，远处回波其实是高空信息，不是地面真实降水。我一般取150km，强对流分析时取100km，只有做拼图评估才用到230km。这个值直接关联画布比例尺，不要用文件里默认的最大量程。

最后是dBZ范围。固定显示范围比自动缩放重要得多。自动缩放会让一场小雨和一个超级单体单胞的色标范围不同，导致颜色与强度对不上，后续人工看图会误判。业务上固定-20到70dBZ，如果某次过程没有那么强，颜色层次少但可比性强。ShowRadarData_intensity里找vmin和vmax并写死，别留成None。

4. 强度图画完不对劲？PPI常见问题与避坑清单

4.1 径向放射状亮线：地物杂波还是解析错位

现象：PPI图上出现从中心向外的一条或几条亮色射线，强度值动辄50dBZ以上，而且固定出现在相同方位角上，不随时间变化。原因多半是地物杂波，雷达波束打到附近高楼、山体或者风力发电机时产生强回波，多普勒雷达在零速度附近的地物特征明显，但强度图上它就是一条线。另一种原因是你解析时读错了径向头，方位角字段没读对，导致某一条径向被重复绘制，也会呈放射状。

解决：先看是不是固定方位。用脚本输出每个方位角的径向平均强度，找到持续超过45dBZ的方位，把它们标记为地物。更严格的办法是使用“径向速度接近0”这个特征，地物回波的径向速度通常在0附近，如果文件同时包含速度数据，就用强度与速度双条件过滤。我的经验是：只在最低仰角做地物剔除，高仰角再做会误删降水。

4.2 大片负值或缺测空洞：字节序和特殊标记的坑

现象：整个PPI有一半区域颜色是-30dBZ以下的蓝色，或者回波内部出现不规则的黑色空洞，看起来像奶酪。

原因：一是字节序读反。信号强度码用两个字节存储时，小端文件用大端读取会让每个值被交换，结果数值分布会在低位和高位之间反复横跳，画出来就是一整片噪声。二是没有处理保留值，某些文件把强度码0当成有效数据参与换算，得到-33dBZ，再经过插值就成了大片蓝色背景。三是插值算法把边缘外的空洞填成低值，看图上远处蓝色区域就知道是投影时网格范围超出数据范围。

解决：先打印原始强度码的直方图，看看值分布是否集中。如果出现大量0和1，就设置掩膜；如果分布均匀且无明显峰值，去改字节序。插值空洞可以用griddata的fill_value=np.nan，然后在绘图时把np.nan设为透明，不要让它参与色标映射。这一条是真正的血泪经验，我见过有人为了避开调字节序，硬把整张图反转颜色，结果把回波强区弄到错误位置。

4.3 强度普遍偏移几dB：换算截距不是别人传的-33

现象：同一场降雨，用脚本画的回波强核位置和业务软件差不多，但整体强度低5dBZ或高5dBZ，颜色对比起来总是差一档。

原因：强度码到dBZ的换算截距不对。网上流传的-33是通用参考值，但雷达厂家在出厂标定时会给每个通道做校正，有的雷达把负增益补偿加在这个截距里，结果截距变成-28或-38。如果你手里的脚本直接写死-33，画出来的强度就会整体平移。

解决：找一份该雷达最近的标定数据，或者对比官方基数据里强度码和业务产品反射率值。取几个已知位置的像素做线性拟合：用文件强度码作为x，业务产品dBZ作为y，拟合斜率k和截距b，然后替换脚本里的参数。没有官方产品时，可以用地面自动站的雨强反推反射率，但只建议验证强降水区，弱回波受衰减影响大。

4.4 回波整体旋转或镜像：方位角坐标系没有对齐

现象：回波结构形态很正常，但是整个图形顺时针旋转了90度，或者左右颠倒，和地图拼图完全对不上。

原因：极坐标投影时，方位角与坐标轴的对应关系搞错了。有的雷达文件里方位角0度指向正北，但绘图公式里x = r * sin(azimuth)、y = r * cos(azimuth)默认0度指向正南，还需要加一个偏北的旋转量；反之，如果文件方位角是顺时针从正北算起，而绘图函数用了逆时针，就会镜像。另一个常见错误是符号问题，正弦和余弦的用法搞反会出现旋转而不是镜像。

解决：先找一个固定地物回波来做方位校正。比如雷达西北方向30公里处有一片明显山体杂波，在图上应该出现在方位角约315度（北北西）的位置。通过旋转角度和镜像开关反复验证，把校正角度写入配置。这种方法比纯看经纬度推算更可靠。ShowRadarData_intensity如果内置方位校正参数，通常在投影函数里有一个azimuth_offset，默认0，你需要在初始化时设成实际体扫基数据的方位角起始值。

5. 批量处理体扫数据：从单张PPI到全套仰角产品

5.1 遍历全部仰角的代码写法与维度陷阱

单仰角跑通后，批量出全部仰角就只是循环的问题，但维度顺序能坑你半天。不同解析库返回的数据维度安排不一样，有的把仰角放在第一维，有的把径向放在第一维。我在写循环前必做一件事：打印数组shape和仰角列表。

import numpy as np from show_intensity import load_radar_file file_path = 'Z_RADR_I_20240815_0200.dat' data, meta = load_radar_file(file_path) print(data.shape) # 期望 (14, 360, 920) 表示 14个仰角 x 360条径向 x 920库 # 如果是 (360, 14, 920)，说明第二维才是仰角

逻辑说明：data.shape[0]是不是仰角数，看仰角列表长度是否等于meta['elevation']的长度。如果维度不对，用np.transpose把仰角轴挪到前面。我吃过一次亏：循环遍历第一维，结果画出来的14张图其实是同一条方位角在不同仰角下的扇面，当时没细看，差点发给预报员。

参数说明：遍历时只输出指定仰角可以省时间。业务上做体扫产品通常选0.5、1.5、2.4、3.4四个仰角，其余仰角用于垂直剖面，不需要每张都出。循环里加一个if elev_idx not in target:跳过，能减少一半以上的IO和画图时间。

5.2 最简单的去地物策略：中值滤波加阈值掩膜

批量处理时不能依赖人工去地物，至少要跑一个最简单的自动策略。中值滤波对孤立地物点效果很好，它的原理是把一个像素替换成周围像素的中值，那些比邻居高20dBZ以上的孤立点会被抹掉。

from scipy.ndimage import median_filter def quality_control(dbz): dbz = np.where(np.isnan(dbz), 0, dbz) filtered = median_filter(dbz, size=3) # 保留强度大于5dBZ且滤波后差异不大的区域 mask = (dbz > 5) & (np.abs(dbz - filtered) < 10) return np.where(mask, dbz, np.nan)

逻辑说明：把缺测值临时填0后做中值滤波，再用原始值与滤波值之差来判断。如果原始值比附近高10dBZ以上，很可能是地物杂波，置为缺测。size=3表示3x3邻域，库长1km时相当于对3km范围做平滑。这个策略简单但有效，缺点是会把小尺度强核也削平。

参数说明：阈值5和10不是固定最优。如果你所在地区地物多，比如山区，可以把5提高到10，把10放宽到15；如果是平原雷达，5和10已经足够。注意只对低仰角做这个操作，高仰角回波本来就是孤立细胞，再做滤波会把真实降水抹掉。

5.3 输出文件命名和解析结果归档

批量出图后，命名规范决定后续检索效率。我用的是气象常规文件名风格：

Z_RADR_I_20240815020000_E001.png Z_RADR_I_20240815020000_E002.png

其中E001对应仰角序号，而不是仰角度数，因为不同雷达仰角列表不一致，用序号最稳。文件名里的时间是雷达开始体扫的时间戳，通常来自文件名本身，不要用图像生成时间，否则没法跟原始基数据对应。

另外，一边画图一边把解析好的三维数组存储为NumPy的.npz格式。格式化的数组比原始基数据小很多，而且后续做ZDR、PHI等产品时不必重新解析。我一般存成：

np.savez_compressed('Z_RADR_I_20240815020000.npz', data=data, meta=meta)

这里的meta是字典，会被自动压成结构化数据。这个习惯帮我躲过好几次原始文件损坏的灾难。基数据可能被覆盖，但解析结果和参数配置留档，就还有后悔药。

6. 少走半年弯路的验证技巧：把强度图叠加地图再跟实况对比

6.1 叠加地图的两种实现

画出来的PPI如果只有色块没有地理参照，很难判断对不对。我最早用的方法是把图存成透明背景PNG，再拖进GIS里跟影像叠加，缺点是不能批量验证。后来改用Python直接叠加地图。

常见做法是用cartopy画海岸线和省界，把雷达PPI的pcolormesh对象放在同一坐标系下。这个方法只需要加一段底图代码，之后所有批量图都可以自动叠加。如果只做内网离线使用，basemap也依然稳定。

import matplotlib.pyplot as plt import cartopy.crs as ccrs import cartopy.feature as cfeature proj = ccrs.PlateCarree() ax = plt.axes(projection=proj) ax.set_extent([118, 123, 29, 34], crs=proj) ax.add_feature(cfeature.COASTLINE, lw=0.6) ax.add_feature(cfeature.BORDERS, linestyle='--', alpha=0.7) # 这里再绘制ppimap = ax.pcolormesh(..., cmap=cmap, norm=norm) plt.savefig('ppi_map_validation.png', dpi=150)

逻辑说明：先确定地图范围，再把雷达中心经纬度换算到地图坐标，PPI数据中的x、y已经是经纬度网格，所以直接传入pcolormesh。cartopy会负责投影转换，不用自己处理。叠加之后重点看两件事：回波形态和地图上的河流、山脉走向是否匹配；强回波中心是否落在城区或特定地形附近。

6.2 用一次真实降水过程给PPI做“体检”

验证脚本是否正确，最有效的办法是找一场你已经知道结果的强降水过程。比如你所在地区发生过一次飑线天气，当时的雷达拼图和雨量计记录还在手里。用ShowRadarData_intensity重跑那段时间的基数据，观察三处：回波带走向是否与拼图一致，强回波中心位置是否与大雨量站重合，不同仰角之间的回波垂直结构是否合理。

如果所有这些都吻合，说明解析、投影、色标三个环节都过关了。只要有一处不对，就回头查对应环节。这个方法比任何单元测试都直接，因为雷达数据的真实物理特征就是最强有力的判据。

我现在每次拿到新的雷达数据格式，依然会先找一场历史强降水过程做验证。当年第一次画PPI方位错了90度，就是靠叠加地图看到回波跑到海边山地才定位到方位角坐标系的问题。这个习惯让无数个调试夜晚少走了弯路。希望帮到你。

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