利用PyCINRAD处理、显示天气雷达基数据

只是简单的安装测试了一下SA/SB/SC格式雷达基数据,不知道什么原因,CD格式的数据读取还存在问题。内容基本参照气象家园上的相关帖子。具体的更多细节可以参照气象家园的帖子,后续我会给出链接。本人不是中尺度研究方向,也是抱着学习的态度来安装测试,也推荐大家业务和科研方面使用PyCINRAD。感谢PyCINRAD的开发者@eeeee@高空急流,花费大量的时间和精力来贡献这款开源库。


注:图片较大显示为压缩后的,点击可下载1下载2原图。

安装

建议git clone https://github.com/CyanideCN/PyCINRAD.git,cd PyCINRAD/, python setup.py install方式来安装,PyCINRAD支持Python3.6+的版本,但是在安装的时候,有时候存在比较慢的情况,可以根据进度来中断,单独来使用pip install packagesconda install packages的方式来安装依赖库,譬如,pyresample和shapefile。

再次,也给出我的PyCINRAD环境,Python环境的导出和导入,一键式安装见我之前写的帖子【一键安装气象常用的Python库(大气科学Python库)】,在此不赘述。除了PyCINRAD的GitHub上提供的安装方法之外,也可以采用我的PyCINRAD.yml文件来解决依赖关系。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
name: PyCINRAD
channels:
- conda-forge
- conda-forge/label/cf201901
- defaults
dependencies:
- arm_pyart=1.9.2=py36_0
- asn1crypto=0.24.0=py36_1003
- bzip2=1.0.6=1
- ca-certificates=2019.3.9=hecc5488_0
- cartopy=0.17.0=py36h929c6f0_1004
- certifi=2019.3.9=py36_0
- cffi=1.12.3=py36hccf1714_0
- cftime=1.0.3.4=py36h7eb728f_0
- chardet=3.0.4=py36_1003
- configobj=5.0.6=py_0
- cryptography=2.6.1=py36h212c5bf_0
- cryptography-vectors=2.3.1=py36_1000
- curl=7.62.0=ha441bb4_0
- cython=0.29.7=py36h6de7cb9_0
- freetype=2.10.0=h24853df_0
- geos=3.7.1=h0a44026_1000
- hdf4=4.2.13=h951d187_2
- hdf5=1.10.4=nompi_h5598ddc_1105
- icu=58.2=h0a44026_1000
- idna=2.8=py36_1000
- jpeg=9c=h1de35cc_1001
- kiwisolver=1.1.0=py36h770b8ee_0
- krb5=1.14.6=0
- libblas=3.8.0=8_openblas
- libcblas=3.8.0=8_openblas
- libcurl=7.62.0=h051b688_0
- libcxx=4.0.1=hcfea43d_1
- libcxxabi=4.0.1=hcfea43d_1
- libedit=3.1.20181209=hb402a30_0
- libffi=3.2.1=h475c297_4
- libgfortran=3.0.1=0
- libiconv=1.15=h01d97ff_1005
- liblapack=3.8.0=8_openblas
- libnetcdf=4.6.2=h2c3f975_1
- libpng=1.6.37=h2573ce8_0
- libssh2=1.8.2=hcdc9a53_2
- libtiff=4.0.10=h79f4b77_1001
- libxml2=2.9.9=hd80cff7_0
- libxslt=1.1.32=h33a18ac_1002
- lxml=4.3.3=py36h08abf6f_0
- matplotlib-base=3.0.2=py36hb2d221d_1
- ncurses=6.1=h0a44026_1
- netcdf4=1.4.2=py36h13743db_0
- olefile=0.46=py_0
- openblas=0.3.6=hd44dcd8_1
- openssl=1.1.1b=h01d97ff_2
- owslib=0.17.1=py_0
- pillow=5.3.0=py36hbddbef0_1000
- pip=19.1.1=py36_0
- proj4=5.2.0=h1de35cc_1001
- pycparser=2.19=py36_1
- pyepsg=0.4.0=py_0
- pykdtree=1.3.1=py36h917ab60_1002
- pyopenssl=19.0.0=py36_0
- pyparsing=2.4.0=py_0
- pyproj=1.9.6=py36h9c430a6_1000
- pyresample=1.12.2=py36h051e8ed_0
- pysocks=1.7.0=py36_0
- python=3.6.8=haf84260_0
- python-dateutil=2.8.0=py_0
- pytz=2019.1=py_0
- pyyaml=5.1=py36h1de35cc_0
- readline=7.0=h1de35cc_5
- requests=2.21.0=py36_1000
- scipy=1.2.1=py36hbd7caa9_1
- setuptools=41.0.1=py36_0
- shapely=1.6.4=py36h4b8df73_1004
- six=1.12.0=py36_1000
- sqlite=3.28.0=ha441bb4_0
- tk=8.6.8=ha441bb4_0
- tornado=6.0.2=py36h01d97ff_0
- trmm_rsl=1.49=3
- urllib3=1.24.2=py36_0
- wheel=0.33.2=py36_0
- xz=5.2.4=h1de35cc_4
- yaml=0.1.7=h1de35cc_1001
- zlib=1.2.11=h1de35cc_3
- pip:
- cinrad==1.3.3
- cycler==0.10.0
- decorator==4.4.0
- dill==0.2.9
- inspyred==1.0.1
- ipython-genutils==0.2.0
- matplotlib==3.1.0rc2
- metpy==0.10.0
- numpy==1.14.5
- packaging==19.0
- pandas==0.24.2
- pint==0.9
- pooch==0.2.1
- pykriging==0.2.0
- pyqt5==5.12.2
- pyqt5-sip==4.19.17
- pyshp==2.1.0
- traitlets==4.3.2
- xarray==0.12.1
- xlrd==1.2.0
prefix: /Users/zhpfu/anaconda3/envs/PyCINRAD

数据、脚本

链接:https://pan.baidu.com/s/1riGd4WBZbXCgggkqct_Eyw 密码:k15x

Z_RADR_I_Z9250_20160701001000_O_DOR_SA_CAP.bin
Z9250_20160701001005_0.6_230_REF_31.2N117.0E_32.4N120.0E.png
Z9250_20160701001005_0.6_230_REF.png

脚本:源于气象家园,版权属于@小其其格@eeeee

  • PPI.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
# -*- coding: utf-8 -*-
"""
Created on Sun Jan 13 21:51:35 2019

@author: Fulang WU
"""

"""定义一个函数,用于读取目标文件夹下的多个雷达基数据"""
#读取指定文件夹中的指定文件类型的文件名
import os

def get_filename(path,fileType):
file_name =[]
final_name = []
for files in os.listdir(path): #root为目录路径 #dirs为路径下的子目录 #files为路径下的所有非子目录
if fileType in files:
file_name.append(files.replace(fileType,''))#生成不带‘.bin’后缀的文件名组成的列表
final_name = [path +item +fileType for item in file_name]#生成‘.bin’后缀的文件名组成的绝对路径列表
return final_name #输出列表
"""-----------------------------------------------------------------------------------------------------------"""


"""读取目标文件夹下的多个雷达基数据"""
Dir = "/Users/zhpfu/Dropbox/Code_Fortress/00_My_Python_Library/Radar_data_process/Method_2/data/" #目标文件夹
fileType = '.bin' #雷达基数据后缀名
sa_radar_file = get_filename(Dir,fileType)
print(sa_radar_file)
"""----------------------------------------------------------------------"""


"""---------------------------绘制雷达图像-------------------------------"""
import cinrad
from cinrad.io import CinradReader, StandardData
from cinrad.visualize import PPI
import matplotlib.ticker as mticker
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER


ele = 0 #第1个仰角
radius = 230 #绘制图像的范围大小

for nFiles in range(len(sa_radar_file)):

f = CinradReader(sa_radar_file[nFiles]) #老版本数据
r = f.get_data(ele,radius,'REF') #0代表第1个仰角
fig = cinrad.visualize.PPI(r) #绘制PPI图像

fig.plot_range_rings(radius, color='white', linewidth=1.0) #绘制圆圈
for i in range(0,radius-30,50):
fig.plot_range_rings(i, color='white', linewidth=1.0) #绘制圆圈

"""
设置经纬度
"""
liner = fig.geoax.gridlines(draw_labels=True,linewidth=2, color='gray', alpha=0.5, linestyle='--')
liner.xlabels_top = False
liner.ylabels_right = False
liner.xformatter = LONGITUDE_FORMATTER
liner.yformatter = LATITUDE_FORMATTER
liner.xlabel_style = {'size': 18, 'color': 'white'}
#liner.xlabel_style = {'color': 'red', 'weight': 'bold'}
liner.ylabel_style = {'size': 18, 'color': 'white'}

#保存图像
fig(Dir)
# /Users/zhpfu/anaconda3/envs/PyCINRAD/lib/python3.6/site-packages/cinrad-1.3.3-py3.6.egg/cinrad/visualize
  • PPI.VCS.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
# -*- coding: utf-8 -*-
"""
Created on Sun Jan 13 21:51:35 2019

@author: Fulang WU
"""

"""定义一个函数,用于读取目标文件夹下的多个雷达基数据"""
#读取指定文件夹中的指定文件类型的文件名
import os
def get_filename(path,fileType):
file_name =[]
final_name = []
for files in os.listdir(path): #root为目录路径 #dirs为路径下的子目录 #files为路径下的所有非子目录
if fileType in files:
file_name.append(files.replace(fileType,''))#生成不带‘.bin’后缀的文件名组成的列表
final_name = [path +item +fileType for item in file_name]#生成‘.bin’后缀的文件名组成的绝对路径列表
return final_name #输出列表
"""-----------------------------------------------------------------------------------------------------------"""



"""读取目标文件夹下的多个雷达基数据"""
Dir = "/Users/zhpfu/Dropbox/Code_Fortress/00_My_Python_Library/Radar_data_process/Method_2/data/" #目标文件夹
fileType = '.bin' #雷达基数据后缀名
#sa_radar_file = Dir+"Z_RADR_I_Z9571_20180304115800_O_DOR_SA_CAP.bin"
sa_radar_file = get_filename(Dir,fileType) #读取
print(sa_radar_file)
"""----------------------------------------------------------------------"""

"""---------------------------绘制雷达图像-------------------------------"""
import cinrad
from cinrad.io import CinradReader, StandardData
from cinrad.visualize import PPI
import matplotlib.ticker as mticker
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER


ele = 0 #第1个仰角
radius = 230 #绘制图像的范围大小

for nFiles in range(len(sa_radar_file)):

f = CinradReader(sa_radar_file[nFiles]) #老版本数据
rl = [f.get_data(i, 230, 'REF') for i in f.angleindex_r]
vcs = cinrad.easycalc.VCS(rl)
sec = vcs.get_section(start_cart=(117.0, 31.2), end_cart=(120., 32.4)) # 传入经纬度坐标
fig = cinrad.visualize.PPI(rl[ele]) #0代表第1个仰角

fig.plot_range_rings(radius, color='white', linewidth=1.0) #绘制圆圈
for i in range(0,radius-30,50):
fig.plot_range_rings(i, color='white', linewidth=1.0) #绘制圆圈

fig.plot_cross_section(sec) #绘制垂直剖面

"""
设置经纬度
"""
liner = fig.geoax.gridlines(draw_labels=True,linewidth=2, color='gray', alpha=0.5, linestyle='--')
liner.xlabels_top = False
liner.ylabels_right = False
liner.xformatter = LONGITUDE_FORMATTER
liner.yformatter = LATITUDE_FORMATTER
liner.xlabel_style = {'size': 18, 'color': 'white'}
#liner.xlabel_style = {'color': 'red', 'weight': 'bold'}
liner.ylabel_style = {'size': 18, 'color': 'white'}

#保存图像
fig(Dir)

绘图

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
(PyCINRAD) ➜  Method_2 python PPI.py
['/Users/zhpfu/Dropbox/Code_Fortress/00_My_Python_Library/Radar_data_process/Method_2/data/Z_RADR_I_Z9250_20160701001000_O_DOR_SA_CAP.bin']
(366, 230)
(366, 230)


(PyCINRAD) ➜ Method_2 PYTHON PPI.VCS.py
['/Users/zhpfu/Dropbox/Code_Fortress/00_My_Python_Library/Radar_data_process/Method_2/data/Z_RADR_I_Z9250_20160701001000_O_DOR_SA_CAP.bin']
(366, 230)
(366, 230)
(367, 230)
(367, 230)
(372, 230)
(372, 230)
(372, 230)
(372, 230)
/Users/zhpfu/anaconda3/envs/PyCINRAD/lib/python3.6/site-packages/cinrad-1.3.3-py3.6-macosx-10.7-x86_64.egg/cinrad/io/io.py:310: UserWarning: Requested data range exceed max range in this tilt
warnings.warn('Requested data range exceed max range in this tilt')
(372, 230)
(372, 230)
(368, 230)
(368, 230)
(365, 230)
(365, 230)
(363, 230)
(363, 230)
(360, 230)
(360, 230)

GUI可视化界面

运行方式:
进到PyCINRAD包的ui文件夹📂:
python main_ui.pyw

打开后,一顿点点点,选择雷达基数据,需要注意的是必须文件名类似于“Z_RADR_I_Z9250_20160701001000_O_DOR_SA_CAP.bin”,不然无法读取。
可视化界面不是很稳定,比较简单的功能!

最后

PyCINRAD其实提供了很丰富的接口,除了可视化之外,对于数据的提取、网格化、分辨率调整都比较灵活,能够支持主要的天气雷达。功能完备!值得深入研究,无论是业务还是科研,都是良好的利器!

问题报错回顾

Mac上不设置是不支持中文输出,有几个关键点需要注意下:
1.在visualize文件下utils.py添加

1
2
3
4
[code=python]import matplotlib.pyplot as plt
plt.rcParams['font.sans-serif'] = ['SimHei'] # 用来正常显示中文标签
plt.rcParams['axes.unicode_minus'] = False # 用来正常显示负号
plt.rcParams['font.family']=['SimHei'][/code]

2.StationNames.xlsx文件确保是包含中文站点名称的(附件名称修改为StationNames.xlsx);并且放到data目录下;
PPI的RAD站点名就可以正常显示了。

3.针对提示:

1
2
3
/Users/zhpfu/anaconda3/envs/PyCINRAD/lib/python3.6/importlib/_bootstrap.py:219: RuntimeWarning: numpy.dtype size changed, may indicate binary incompatibility. Expected 96, got 88
return f(*args, **kwds)
Cython is not installed, velocity dealias function cannot be used

解决办法:
pip uninstall cinrad
再从GitHub上git clone最新的文件,进到文件夹内:
python setup.py install
即可。

4.针对警告:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
/Users/zhpfu/anaconda3/envs/PyCINRAD/lib/python3.6/importlib/_bootstrap.py:219: RuntimeWarning: numpy.dtype size changed, may indicate binary incompatibility. Expected 96, got 88
return f(*args, **kwds)

解决办法很简单,使用更低版本numpy即可:
(PyCINRAD) ➜ pip install numpy==1.14.5
Looking in indexes: https://pypi.tuna.tsinghua.edu.cn/simple
Collecting numpy==1.14.5
Downloading https://pypi.tuna.tsinghua.edu.cn/packages/f6/cd/b2c50b5190b66c711c23ef23c41d450297eb5a54d2033f8dcb3b8b13ac85/numpy-1.14.5-cp36-cp36m-macosx_10_6_intel.macosx_10_9_intel.macosx_10_9_x86_64.macosx_10_10_intel.macosx_10_10_x86_64.whl (4.7MB)
|████████████████████████████████| 4.7MB 714kB/s
Installing collected packages: numpy
Found existing installation: numpy 1.15.0
Uninstalling numpy-1.15.0:
Successfully uninstalled numpy-1.15.0
Successfully installed numpy-1.14.5

以上是一些问题解决办法!仅供参考!

问题:PyCINRAD目前还不支持CD的PPI绘制?报错如下;

1
2
3
4
5
6
7
8
9
10
Traceback (most recent call last):
File "PPI_CD.py", line 36, in <module>
r = f.get_data(5, 100, 'REF')
File "/Users/zhpfu/anaconda3/envs/PyCINRAD/lib/python3.6/site-packages/cinrad-1.3.3-py3.6.egg/cinrad/io/io.py", line 360, in get_data
x, y, z, d, a = self.projection(reso)
File "/Users/zhpfu/anaconda3/envs/PyCINRAD/lib/python3.6/site-packages/cinrad-1.3.3-py3.6.egg/cinrad/io/io.py", line 371, in projection
lonx, latx = get_coordinate(r, theta, self.elev, self.stationlon, self.stationlat, h_offset=h_offset)
File "/Users/zhpfu/anaconda3/envs/PyCINRAD/lib/python3.6/site-packages/cinrad-1.3.3-py3.6.egg/cinrad/projection.py", line 62, in get_coordinate
deltav = np.cos(azimuth) * distance * np.cos(elev * deg2rad)
AttributeError: 'NoneType' object has no attribute 'cos'

参考

1.网址1
2.网址2
3.网址3
4.网址4

有任何问题都欢迎交流探讨,共同学习进步!


-------------本文结束 感谢阅读-------------
作者Gavin
有问题请在相应页面留言(评论系统DISQUS需要"翻墙"才可见)
或者私信给我 GMAIL: zhpfu.atm@gmail.com
满分是10分的话,这篇文章你给几分
--> -->