Python可视化 | 中尺度对流系统反射率截面

2021-08-26 18:21:54 浏览数 (1)

背景介绍

中尺度对流系统,简称MCS(Mesoscale Convective System),是造成暴雨 、冰雹 、雷雨大风和龙卷等灾害性天气的重要系统。

Orlanski(1975)按尺度将MCS划分为α、β、γ三种中尺度,α中尺度对流系统(MαCS)水平尺度为200~2000km,β中尺度对流系统(MβCS)为20~200km,γ中尺度对流系统(MγCS)为2~20km。

按对流系统的组织形式分为三类:孤立对流系统、带状对流系统以及圆形对流系统或中尺度对流复合体MCC(Mesoscale Convective Complex)。

孤立对流系统有三种类型:⑴普通单体风暴;⑵多单体风暴;⑶超级单体风暴。带状对流系统最典型的代表就是飑线系统。

导入模块

代码语言:javascript复制
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import matplotlib as mpl
import matplotlib.colors as colors
import matplotlib.pyplot as plt
import numpy as np
import xarray as xr
import metpy.calc as mpcalc
from metpy.units import units
from metpy.interpolate import cross_section, log_interpolate_1d

import warnings
warnings.filterwarnings('ignore')

自定义色阶

代码语言:javascript复制
def yuv_rainbow_24(nc):
    path1 = np.linspace(0.8*np.pi, 1.8*np.pi, nc)
    path2 = np.linspace(-0.33*np.pi, 0.33*np.pi, nc)

    y = np.concatenate([np.linspace(0.3, 0.85, nc*2//5),
                        np.linspace(0.9, 0.0, nc - nc*2//5)])
    u = 0.40*np.sin(path1)
    v = 0.55*np.sin(path2)   0.1

    rgb_from_yuv = np.array([[1, 0, 1.13983],
                             [1, -0.39465, -0.58060],
                             [1, 2.03211, 0]])
    cmap_dict = {'blue': [], 'green': [], 'red': []}
    for i in range(len(y)):
        yuv = np.array([y[i], u[i], v[i]])
        rgb = rgb_from_yuv.dot(yuv)
        red_tuple = (i/(len(y)-1), rgb[0], rgb[0])
        green_tuple = (i/(len(y)-1), rgb[1], rgb[1])
        blue_tuple = (i/(len(y)-1), rgb[2], rgb[2])
        cmap_dict['blue'].append(blue_tuple)
        cmap_dict['red'].append(red_tuple)
        cmap_dict['green'].append(green_tuple)

    return cmap_dict

cmap = colors.LinearSegmentedColormap('homeyer_rainbow', yuv_rainbow_24(15), mpl.rcParams['image.lut'])

读取数据

代码语言:javascript复制
data = xr.open_dataset('/home/mw/input/data3672/wrf_for_metpy_demo.nc').metpy.parse_cf()
data

反射率绘制

代码语言:javascript复制
data['dbz'].max(data['dbz'].metpy.vertical.name).plot(cmap=cmap, vmin=0,
                                                      vmax=60, figsize=(14,6.),
                                                      aspect='equal')
代码语言:javascript复制
<matplotlib.collections.QuadMesh at 0x7f2c15ca67d0>

确定端点位置

代码语言:javascript复制
start, end = ((48.423, -99.771),(43.476, -88.490))

获取截面数据

代码语言:javascript复制
cross = cross_section(data, start, end, 1000)
cross

绘制中尺度对流截面

代码语言:javascript复制
fig = plt.figure(1, figsize=(14., 6.))
ax = plt.axes()
pres_axis = np.arange(1000, 80, -10) * units.hPa

with warnings.catch_warnings():
    warnings.simplefilter("ignore")
    dbz_interp = log_interpolate_1d(pres_axis, cross['pressure'], cross['dbz'].values)

refl_contour = ax.contourf(cross['XLONG'], pres_axis, dbz_interp,
                           levels=np.arange(0, 75, 5), cmap=cmap)
refl_colorbar = fig.colorbar(refl_contour)

ax.set_yscale('symlog')
ax.set_yticklabels(np.arange(1000, 50, -100))
ax.set_ylim(pres_axis.max(), pres_axis.min())
ax.set_yticks(np.arange(1000, 50, -100))

data_crs = data['dbz'].metpy.cartopy_crs
ax_inset = fig.add_axes([0.072, 0.63, 0.25, 0.25], projection=data_crs)


ax_inset.contourf(data['west_east'], data['south_north'], data['dbz'].max('bottom_top'),
                  levels=np.arange(0, 75, 5), cmap=cmap, zorder=5)


endpoints = data_crs.transform_points(ccrs.Geodetic(),
                                      *np.vstack([start, end]).transpose()[::-1])
ax_inset.scatter(endpoints[:, 0], endpoints[:, 1], c='k', zorder=10)
ax_inset.plot(cross['west_east'], cross['south_north'], c='k', zorder=10)

ax_inset.add_feature(cfeature.LAND.with_scale('50m'), facecolor=cfeature.COLORS['land'], zorder=1)
ax_inset.add_feature(cfeature.STATES.with_scale('50m'), edgecolor='#9f9f68', zorder=2)
ax_inset.add_feature(cfeature.LAKES.with_scale('50m'), facecolor=cfeature.COLORS['water'],
               edgecolor='#9f9f68', zorder=2)

ax_inset.set_title('')
ax.set_title('WRF Simulated Reflectivity - 2015-07-13 05:00Z Forecastn'
             'Initialized 2015-07-12 12:00Z, Thompson Microphysics, MYJ PBLn'
             'Inset: Composite Reflectivity')
ax.set_ylabel('Pressure (hPa)')
ax.set_xlabel('Longitude (degrees east)')
refl_colorbar.set_label('Z (dbZ)')

plt.show()

获取数据

代码及数据获取请在好奇心Log公众号后台回复wrf_demo,数据1.7GB,可能耗时较长

更多WRF示例数据:https://www.heywhale.com/mw/dataset/606575296a76ca0017484b60

0 人点赞