两种方法将气象雷达数据转为易处理的格式
个人信息
公众号：气python风雨
关注我获取更多学习资料，第一时间收到我的Python学习资料，也可获取我的联系方式沟通合作
温馨提示
由于可视化代码过长隐藏，可点击运行Fork查看
若没有成功加载可视化图，点击运行可以查看
ps:隐藏代码在【代码已被隐藏】所在行，点击所在行，可以看到该行的最右角，会出现个三角形，点击查看即可
前言
项目目标
读者来信：我想获得一个雷达基数据的里每个有效数据点的反射率强度，经纬度，海拔高度，这样一个三维的反射率强度数据，我想找到反射率强度达到某个值的这个或者这一组点从中心最强到临近区域最弱区域的三维距离和梯度变化，主要是想看看山区局地的小对流云团的特征。
项目方法
在以下内容中，展示两种方法分别将雷达数据转为易于处理的表格数据和三维xarray数据
!pip install --upgrade --index-url=https://pypi.mirrors.ustc.edu.cn/simple cinrad --userLooking in indexes: https://pypi.mirrors.ustc.edu.cn/simple
Requirement already satisfied: cinrad in /opt/conda/lib/python3.9/site-packages (1.8.0)
Collecting cinrad
Downloading https://mirrors.ustc.edu.cn/pypi/packages/5a/5e/6144ebcccd0e65849b533e931de6ff894c2c4c4854cc667705da3e2198ef/cinrad-1.9.1.tar.gz (395 kB)
[2K [90m━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━[0m [32m395.0/395.0 kB[0m [31m5.8 MB/s[0m eta [36m0:00:00[0ma [36m0:00:01[0m
[?25h Preparing metadata (setup.py) ... [?25ldone
[?25hRequirement already satisfied: metpy>=0.8 in /opt/conda/lib/python3.9/site-packages (from cinrad) (1.6.1)
Requirement already satisfied: cartopy>=0.15 in /opt/conda/lib/python3.9/site-packages (from cinrad) (0.23.0)
Requirement already satisfied: pyshp!=2.0.0,!=2.0.1 in /opt/conda/lib/python3.9/site-packages (from cinrad) (2.3.1)
Requirement already satisfied: matplotlib>=2.2 in /opt/conda/lib/python3.9/site-packages (from cinrad) (3.8.3)
Requirement already satisfied: vanadis in /opt/conda/lib/python3.9/site-packages (from cinrad) (0.0.2)
Requirement already satisfied: cinrad_data>=0.1 in /opt/conda/lib/python3.9/site-packages (from cinrad) (0.1)
Requirement already satisfied: shapely>=1.7 in /opt/conda/lib/python3.9/site-packages (from cartopy>=0.15->cinrad) (1.8.5.post1)
Requirement already satisfied: packaging>=20 in /opt/conda/lib/python3.9/site-packages (from cartopy>=0.15->cinrad) (23.2)
Requirement already satisfied: numpy>=1.21 in /opt/conda/lib/python3.9/site-packages (from cartopy>=0.15->cinrad) (1.26.4)
Requirement already satisfied: pyproj>=3.3.1 in /opt/conda/lib/python3.9/site-packages (from cartopy>=0.15->cinrad) (3.4.1)
Requirement already satisfied: contourpy>=1.0.1 in /opt/conda/lib/python3.9/site-packages (from matplotlib>=2.2->cinrad) (1.2.0)
Requirement already satisfied: pillow>=8 in /opt/conda/lib/python3.9/site-packages (from matplotlib>=2.2->cinrad) (9.4.0)
Requirement already satisfied: pyparsing>=2.3.1 in /opt/conda/lib/python3.9/site-packages (from matplotlib>=2.2->cinrad) (3.0.9)
Requirement already satisfied: cycler>=0.10 in /opt/conda/lib/python3.9/site-packages (from matplotlib>=2.2->cinrad) (0.11.0)
Requirement already satisfied: python-dateutil>=2.7 in /opt/conda/lib/python3.9/site-packages (from matplotlib>=2.2->cinrad) (2.8.2)
Requirement already satisfied: importlib-resources>=3.2.0 in /opt/conda/lib/python3.9/site-packages (from matplotlib>=2.2->cinrad) (5.7.1)
Requirement already satisfied: fonttools>=4.22.0 in /opt/conda/lib/python3.9/site-packages (from matplotlib>=2.2->cinrad) (4.33.3)
Requirement already satisfied: kiwisolver>=1.3.1 in /opt/conda/lib/python3.9/site-packages (from matplotlib>=2.2->cinrad) (1.4.2)
Requirement already satisfied: xarray>=0.21.0 in /opt/conda/lib/python3.9/site-packages (from metpy>=0.8->cinrad) (2024.2.0)
Requirement already satisfied: pandas>=1.4.0 in /opt/conda/lib/python3.9/site-packages (from metpy>=0.8->cinrad) (2.0.3)
Requirement already satisfied: scipy>=1.8.0 in /opt/conda/lib/python3.9/site-packages (from metpy>=0.8->cinrad) (1.11.4)
Requirement already satisfied: traitlets>=5.0.5 in /opt/conda/lib/python3.9/site-packages (from metpy>=0.8->cinrad) (5.2.0)
Requirement already satisfied: pooch>=1.2.0 in /opt/conda/lib/python3.9/site-packages (from metpy>=0.8->cinrad) (1.8.1)
Requirement already satisfied: pint>=0.17 in /opt/conda/lib/python3.9/site-packages (from metpy>=0.8->cinrad) (0.23)
Requirement already satisfied: zipp>=3.1.0 in /opt/conda/lib/python3.9/site-packages (from importlib-resources>=3.2.0->matplotlib>=2.2->cinrad) (3.8.0)
Requirement already satisfied: pytz>=2020.1 in /opt/conda/lib/python3.9/site-packages (from pandas>=1.4.0->metpy>=0.8->cinrad) (2022.1)
Requirement already satisfied: tzdata>=2022.1 in /opt/conda/lib/python3.9/site-packages (from pandas>=1.4.0->metpy>=0.8->cinrad) (2024.1)
Requirement already satisfied: typing-extensions in /opt/conda/lib/python3.9/site-packages (from pint>=0.17->metpy>=0.8->cinrad) (4.7.1)
Requirement already satisfied: requests>=2.19.0 in /opt/conda/lib/python3.9/site-packages (from pooch>=1.2.0->metpy>=0.8->cinrad) (2.27.1)
Requirement already satisfied: platformdirs>=2.5.0 in /opt/conda/lib/python3.9/site-packages (from pooch>=1.2.0->metpy>=0.8->cinrad) (4.2.0)
Requirement already satisfied: certifi in /opt/conda/lib/python3.9/site-packages (from pyproj>=3.3.1->cartopy>=0.15->cinrad) (2024.2.2)
Requirement already satisfied: six>=1.5 in /opt/conda/lib/python3.9/site-packages (from python-dateutil>=2.7->matplotlib>=2.2->cinrad) (1.16.0)
Requirement already satisfied: charset-normalizer~=2.0.0 in /opt/conda/lib/python3.9/site-packages (from requests>=2.19.0->pooch>=1.2.0->metpy>=0.8->cinrad) (2.0.12)
Requirement already satisfied: idna<4,>=2.5 in /opt/conda/lib/python3.9/site-packages (from requests>=2.19.0->pooch>=1.2.0->metpy>=0.8->cinrad) (3.3)
Requirement already satisfied: urllib3<1.27,>=1.21.1 in /opt/conda/lib/python3.9/site-packages (from requests>=2.19.0->pooch>=1.2.0->metpy>=0.8->cinrad) (1.26.9)
Building wheels for collected packages: cinrad
Building wheel for cinrad (setup.py) ... [?25ldone
[?25h Created wheel for cinrad: filename=cinrad-1.9.1-cp39-cp39-linux_x86_64.whl size=600590 sha256=54b6f6bb4c789ea96849bb19111d915862080ca487326379154519e6c8eba61a
Stored in directory: /home/mw/.cache/pip/wheels/60/8f/52/2f687389af413cd0d5248bb3580cfbafdf5d969502ec7808c8
Successfully built cinrad
Installing collected packages: cinrad
Successfully installed cinrad-1.9.1安装完毕后重启一下内核
方法一：雷达数据转csv
导入pycinrad库与读取
import cinrad
f = cinrad.io.read_auto('/home/mw/input/pycwr5461/Z_RADR_I_Z9898_20190828181529_O_DOR_SAD_CAP_FMT (1).bin.bz2')
f.available_product(0)['TREF', 'REF', 'SQI', 'ZDR', 'RHO', 'PHI', 'KDP', 'SNRH']取第一个仰角反射率数据
ele = 0 #选择第1个仰角
radius = 400 #绘制图像的范围大小，单位km
r = f.get_data(ele,radius,"REF") #选择反射率数据
r<xarray.Dataset> Size: 19MB
Dimensions: (azimuth: 366, distance: 1600)
Coordinates:
• azimuth (azimuth) float32 1kB 4.689 4.707 4.724 ... 4.658 4.675 4.692
• distance (distance) float64 13kB 0.25 0.5 0.75 1.0 ... 399.5 399.8 400.0
Data variables:
REF (azimuth, distance) float64 5MB nan nan nan nan ... nan nan nan
longitude (azimuth, distance) float64 5MB 110.2 110.2 110.2 ... 106.4 106.4
latitude (azimuth, distance) float64 5MB 20.0 20.0 20.0 ... 19.92 19.92
height (azimuth, distance) float64 5MB 0.1201 0.1222 ... 12.89 12.9
Attributes:
elevation: 0.48339844
range: 400
scan_time: 2019-08-28 18:15:29
site_code: Z9898
site_name: 海口
site_longitude: 110.245834
site_latitude: 19.99639
tangential_reso: 0.25
nyquist_vel: 8.83834
task: VCP21D
以上为xarray数据，那么我们很方便将其转为pandas的格式
# 将xarray DataArray转换为pandas DataFrame
df = r.to_dataframe()
print(df) REF longitude latitude height
azimuth distance
4.689176 0.25 NaN 110.243438 19.996337 0.120113
0.50 NaN 110.241042 19.996285 0.122233
0.75 NaN 110.238646 19.996233 0.124361
1.00 NaN 110.236250 19.996180 0.126496
1.25 3.5 110.233854 19.996128 0.128638
... ... ... ... ...
4.692318 399.00 NaN 106.423289 19.924250 12.849046
399.25 NaN 106.420895 19.924204 12.862894
399.50 NaN 106.418501 19.924159 12.876749
399.75 NaN 106.416107 19.924114 12.890612
400.00 NaN 106.413713 19.924069 12.904482
[585600 rows x 4 columns]取全部仰角反射率数据
rl = [f.get_data(i, 400, 'REF') for i in f.angleindex_r]
rl[1]<xarray.Dataset> Size: 19MB
Dimensions: (azimuth: 366, distance: 1600)
Coordinates:
• azimuth (azimuth) float32 1kB 5.646 5.664 5.681 ... 5.614 5.631 5.649
• distance (distance) float64 13kB 0.25 0.5 0.75 1.0 ... 399.5 399.8 400.0
Data variables:
REF (azimuth, distance) float64 5MB nan nan nan nan ... nan nan nan
longitude (azimuth, distance) float64 5MB 110.2 110.2 110.2 ... 107.9 107.9
latitude (azimuth, distance) float64 5MB 20.0 20.0 20.0 ... 22.9 22.9
height (azimuth, distance) float64 5MB 0.1245 0.1311 ... 19.94 19.96
Attributes:
elevation: 1.4941406
range: 400
scan_time: 2019-08-28 18:15:29
site_code: Z9898
site_name: 海口
site_longitude: 110.245834
site_latitude: 19.99639
tangential_reso: 0.25
nyquist_vel: 8.83834
task: VCP21D
以上则是多个仰角xarray数据的列表
import pandas as pd
# 将每个xarray DataArray转换为pandas DataFrame
df_list = [da.to_dataframe() for da in rl]
# 使用pandas concat函数拼接所有的DataFrame
df_concatenated = pd.concat(df_list, ignore_index=True)
print(df_concatenated) REF longitude latitude height
0 NaN 110.243438 19.996337 0.120113
1 NaN 110.241042 19.996285 0.122233
2 NaN 110.238646 19.996233 0.124361
3 NaN 110.236250 19.996180 0.126496
4 3.5 110.233854 19.996128 0.128638
... ... ... ... ...
5236795 NaN 112.333156 22.784978 142.748626
5236796 NaN 112.334491 22.786726 142.843865
5236797 NaN 112.335825 22.788473 142.939111
5236798 NaN 112.337160 22.790220 143.034364
5236799 NaN 112.338495 22.791967 143.129625
[5236800 rows x 4 columns]方法二：使用pyart库进行格点化
pycinrad转pyart格式
from cinrad.io.export import standard_data_to_pyart
radar = standard_data_to_pyart(f)## You are using the Python ARM Radar Toolkit (Py-ART), an open source
## library for working with weather radar data. Py-ART is partly
## supported by the U.S. Department of Energy as part of the Atmospheric
## Radiation Measurement (ARM) Climate Research Facility, an Office of
## Science user facility.
##
## If you use this software to prepare a publication, please cite:
##
## JJ Helmus and SM Collis, JORS 2016, doi: 10.5334/jors.119
/home/mw/.local/lib/python3.9/site-packages/cinrad/io/level2.py:780: RuntimeWarning: Empty data
warnings.warn("Empty data", RuntimeWarning)插值为三维数据
import pyart
grid = pyart.map.grid_from_radars(
(radar,),
grid_shape=(20, 500, 500), # 设置了网格的形状为 (200, 500, 500)，想计算速度快点可以减少格点数
grid_limits=( ( 0.0, 20000,),(-200000.0, 200000.0),(-200000, 200000.0)), #设置了网格的范围
grid_origin=(32.2 ,118.7), # 绘图中心
fields=["reflectivity"], # 变量
)ds = grid.to_xarray()
ds<xarray.Dataset> Size: 64MB
Dimensions: (time: 1, z: 20, y: 500, x: 500)
Coordinates:
• time (time) object 8B 2019-08-28 18:15:29.410594
• z (z) float64 160B 0.0 1.053e+03 2.105e+03 ... 1.895e+04 2e+04
lat (y, x) float64 2MB 30.38 30.38 30.38 ... 33.98 33.98 33.98
lon (y, x) float64 2MB 116.6 116.6 116.6 ... 120.9 120.9 120.9• y (y) float64 4kB -2e+05 -1.992e+05 ... 1.992e+05 2e+05
• x (x) float64 4kB -2e+05 -1.992e+05 ... 1.992e+05 2e+05
Data variables:
reflectivity (time, z, y, x) float64 40MB nan nan nan nan ... nan nan nan
ROI (time, z, y, x) float32 20MB 3.457e+04 3.458e+04 ... 5.01e+04
可以看到数据已经是有xyz的三维数据了，那不是随意拿捏
剩余的计算就自行解决吧
小结
为了实现上述目标，项目采用了两种不同的方法来转换原始雷达数据，使其更便于后续的数据处理与分析：
表格数据转换：首先将雷达基数据转化为表格形式，这样可以方便地使用传统的数据分析工具进行处理。表格数据结构清晰，便于观察单个数据点的各项属性，比如反射率强度、地理位置坐标（经纬度）以及海拔高度等。
三维xarray数据转换：此外，还利用了xarray库将雷达数据组织成三维数据集。xarray是一个Python库，它提供了带有标签的多维数组，非常适合于气象和地理空间数据的存储和操作。通过这种方式，可以在空间维度上（如纬度、经度、高度）直接进行高效的数值计算和数据分析，特别适合于研究反射率强度的空间分布及其梯度变化。
这两种方法各有优势，表格数据更适合直观查看和基础统计分析，而xarray则更适合复杂的多维数据分析和科学计算。通过结合使用这两种方式，可以全面深入地了解雷达数据中的信息
