Introduction to the thundeR package
Bartosz Czernecki, Mateusz Taszarek, Piotr Szuster
Source:vignettes/getstarted.Rmd
getstarted.RmdthundeR
thundeR is a freeware R package and
collection of functions for rapid computation and visualisation of
convective parameters commonly used in the operational forecasting of
severe convective storms. Core algorithm is based on C++ code
implemented into R language via Rcpp. This solution allows
to compute over 200 thermodynamic and kinematic parameters in less than
0.02s per profile and process large datasets such as reanalyses or
operational NWP models in a reasonable amount of time. Package has been
developed since 2017 by research meteorologists specializing in severe
convective storms and is constantly updated with new features.
Online browser:
Online rawinsonde browser of thundeR package is
available at http://rawinsonde.com
Main functions:
sounding_compute()- A core function for calculating convective parameters commonly used in the operational prediction of severe convective storms. Returns a vector of parameters.sounding_plot()- Function to plot a composite of Skew-T, hodograph and selected convective parameters on a single layoutsounding_save()- Auxiliary function tosounding_plot` that plots a composite of Skew-T, hodograph and selected convective parameters on a single layout and exports graphical file.get_sounding()- Download rawinsonde measurement from sounding database of the University of Wyoming in a form convenient to use withthundeRpackage.
Examples
Examples show aplication of selected thundeR package
functions
Example 1
Draw Skew-T, hodograph and convective parameters on a single layout and export to png file:
library(thunder)
data("sounding_vienna") # load example dataset (Vienna rawinsonde profile for 23 Aug 2011 12UTC):
pressure = sounding_vienna$pressure # vector of pressure [hPa]
altitude = sounding_vienna$altitude # vector of altitude [meters]
temp = sounding_vienna$temp # vector of temperature [degree Celsius]
dpt = sounding_vienna$dpt # vector of dew point temperature [degree Celsius]
wd = sounding_vienna$wd # vector of wind direction [azimuth in degrees]
ws = sounding_vienna$ws # vector of wind speed [knots]
sounding_save(filename = "Vienna.png",
title = "Vienna - 23 August 2011 1200 UTC",
pressure, altitude, temp, dpt, wd, ws)
#> Your display device is 10 x 6 in.
#> It is recommended to use at least 10 x 7.5 in. plotting window
#> or consider saving the layout into file
Example 2
Download North Platte rawinsonde profile for 03 Jul 1999 00UTC and export to png file:
profile = get_sounding(wmo_id = 72562, yy = 1999, mm = 7, dd = 3, hh = 0)
head(profile) # show first few rows of downloaded dataset
#> pressure altitude temp dpt wd ws
#> 2 972.0 287 24.0 13.0 130 12
#> 3 970.0 305 24.1 12.8 130 13
#> 4 964.0 359 24.2 12.2 134 16
#> 5 936.6 610 22.2 11.6 150 31
#> 6 925.0 718 21.4 11.4 150 32
#> 7 904.1 914 19.9 10.2 160 33
sounding_save(filename = "NorthPlatte.png", title = "North Platte - 03 July 1999 0000 UTC", profile$pressure, profile$altitude, profile$temp, profile$dpt, profile$wd, profile$ws)
#> Your display device is 10 x 6 in.
#> It is recommended to use at least 10 x 7.5 in. plotting window
#> or consider saving the layout into file
Example 3
Compute convective parameters based on a sample vertical profile data:
pressure = c(1000, 855, 700, 500, 300, 100, 10) # pressure [hPa]
altitude = c(0, 1500, 2500, 6000, 8500, 12000, 25000) # altitude [meters]
temp = c(25, 10, 0, -15, -30, -50, -92) # temperature [degree Celsius]
dpt = c(20, 5, -5, -30, -55, -80, -99) # dew point temperature [degree Celsius]
wd = c(0, 90, 135, 180, 270, 350, 0) # wind direction [azimuth in degress]
ws = c(5, 10, 20, 30, 40, 5, 0) # wind speed [knots]
accuracy = 2 # accuracy of computations where 3 = high (slow), 2 = medium (recommended), 1 = low (fast)
sounding_compute(pressure, altitude, temp, dpt, wd, ws, accuracy)
#> MU_CAPE MU_CAPE_M10 MU_CAPE_M10_PT
#> 2269.925733301 998.144308627 313.074714966
#> MU_02km_CAPE MU_03km_CAPE MU_HGL_CAPE
#> 247.979380820 575.629305978 1616.538411895
#> MU_CIN MU_LCL_HGT MU_LFC_HGT
#> 0.000000000 730.000000000 730.000000000
#> MU_EL_HGT MU_LI MU_LI_M10
#> 8300.000000000 -10.111930240 -10.853947310
#> MU_WMAX MU_EL_TEMP MU_LCL_TEMP
#> 67.378419888 -28.800000000 17.700000000
#> MU_LFC_TEMP MU_MIXR MU_CAPE_500
#> 17.700000000 14.875917910 1076.032178138
#> MU_CAPE_500_M10 MU_CAPE_500_M10_PT MU_CIN_500
#> 389.365140754 137.081449189 0.000000000
#> MU_LI_500 MU_LI_500_M10 SB_CAPE
#> -5.041702392 -6.234631722 2269.925733301
#> SB_CAPE_M10 SB_CAPE_M10_PT SB_02km_CAPE
#> 998.144308627 313.074714966 247.979380820
#> SB_03km_CAPE SB_HGL_CAPE SB_CIN
#> 575.629305978 1616.538411895 0.000000000
#> SB_LCL_HGT SB_LFC_HGT SB_EL_HGT
#> 730.000000000 730.000000000 8300.000000000
#> SB_LI SB_LI_M10 SB_WMAX
#> -10.111930240 -10.853947310 67.378419888
#> SB_EL_TEMP SB_LCL_TEMP SB_LFC_TEMP
#> -28.800000000 17.700000000 17.700000000
#> SB_MIXR ML_CAPE ML_CAPE_M10
#> 14.875917910 1646.063927273 670.100095437
#> ML_CAPE_M10_PT ML_02km_CAPE ML_03km_CAPE
#> 225.281641180 164.079828657 422.429000173
#> ML_HGL_CAPE ML_CIN ML_LCL_HGT
#> 1250.022127396 0.000000000 975.000000000
#> ML_LFC_HGT ML_EL_HGT ML_LI
#> 975.000000000 7900.000000000 -7.620329500
#> ML_LI_M10 ML_WMAX ML_EL_TEMP
#> -8.584526572 57.377067323 -26.400000000
#> ML_LCL_TEMP ML_LFC_TEMP ML_MIXR
#> 15.250000000 15.250000000 13.048747763
#> LR_0500m LR_01km LR_02km
#> -10.000000000 -10.000000000 -10.000000000
#> LR_03km LR_04km LR_06km
#> -9.047619048 -7.857142857 -6.666666667
#> LR_16km LR_26km LR_24km
#> -6.000000000 -5.000000000 -5.767195767
#> LR_36km LR_26km_MAX LR_500700hPa
#> -4.285714286 -5.714285714 -4.285714286
#> LR_500800hPa LR_600800hPa FRZG_HGT
#> -5.180722892 -5.833333333 2500.000000000
#> FRZG_wetbulb_HGT HGT_max_thetae_03km HGT_min_thetae_04km
#> 2275.000000000 0.000000000 3750.000000000
#> Delta_thetae Delta_thetae_min04km Thetae_01km
#> 28.069785504 28.834554624 330.532272184
#> Thetae_02km DCAPE Cold_Pool_Strength
#> 323.619051941 598.310049405 12.632184321
#> Wind_Index PRCP_WATER Moisture_Flux_02km
#> 33.906407176 27.104631918 30.425480722
#> RH_01km RH_02km RH_14km
#> 0.729103473 0.719706866 0.645204350
#> RH_25km RH_36km RH_HGL
#> 0.555017886 0.443603638 0.460317694
#> BS_0500m BS_01km BS_02km
#> 1.917219591 3.834439182 8.782108409
#> BS_03km BS_06km BS_08km
#> 12.655968060 18.005540000 17.407678201
#> BS_36km BS_26km BS_16km
#> 9.369289559 13.330420492 16.647826116
#> BS_18km BS_EFF_MU BS_EFF_SB
#> 20.279072460 14.223187917 14.223187917
#> BS_EFF_ML BS_SFC_to_M10 BS_1km_to_M10
#> 13.896750821 15.510376997 13.649859352
#> BS_2km_to_M10 BS_MU_LFC_to_M10 BS_SB_LFC_to_M10
#> 9.883032095 14.073701016 14.073701016
#> BS_ML_LFC_to_M10 BS_MW02_to_SM BS_MW02_to_RM
#> 13.686365658 7.304027597 10.140964728
#> BS_MW02_to_LM BS_HGL_to_SM BS_HGL_to_RM
#> 10.786958453 4.893411960 7.785950290
#> BS_HGL_to_LM MW_0500m MW_01km
#> 9.988490342 2.308638103 2.425112273
#> MW_02km MW_03km MW_06km
#> 3.347646594 4.800331697 7.810700539
#> MW_13km SRH_100m_RM SRH_250m_RM
#> 6.838886565 4.253473515 10.053664673
#> SRH_500m_RM SRH_1km_RM SRH_3km_RM
#> 19.720649935 39.634639574 152.521946573
#> SRH_36km_RM SRH_100m_LM SRH_250m_LM
#> 236.590130791 1.502688386 3.551808913
#> SRH_500m_LM SRH_1km_LM SRH_3km_LM
#> 6.967009791 14.002323601 -13.130823085
#> SRH_36km_LM SV_500m_RM SV_01km_RM
#> -24.378976558 0.003904164 0.003879744
#> SV_03km_RM SV_500m_LM SV_01km_LM
#> 0.004796086 0.001013916 0.001100521
#> SV_03km_LM MW_SR_500m_RM MW_SR_01km_RM
#> -0.001435868 10.086252358 10.150075094
#> MW_SR_03km_RM MW_SR_500m_LM MW_SR_01km_LM
#> 9.535939066 13.782126361 12.858486067
#> MW_SR_03km_LM MW_SR_VM_500m_RM MW_SR_VM_01km_RM
#> 8.457941382 10.107849533 10.225339777
#> MW_SR_VM_03km_RM MW_SR_VM_500m_LM MW_SR_VM_01km_LM
#> 10.535820052 13.764654240 12.837090856
#> MW_SR_VM_03km_LM SV_FRA_500m_RM SV_FRA_01km_RM
#> 8.834244774 0.998219354 0.987136746
#> SV_FRA_03km_RM SV_FRA_500m_LM SV_FRA_01km_LM
#> 0.956008752 0.259238846 0.280009490
#> SV_FRA_03km_LM Bunkers_RM_A Bunkers_RM_M
#> -0.286213018 209.404566875 7.793284819
#> Bunkers_LM_A Bunkers_LM_M Bunkers_MW_A
#> 122.058537262 13.182518634 151.949436108
#> Bunkers_MW_M Corfidi_downwind_A Corfidi_downwind_M
#> 7.810700539 218.695496108 14.698235841
#> Corfidi_upwind_A Corfidi_upwind_M K_Index
#> 231.328321825 9.179441145 24.354838710
#> Showalter_Index TotalTotals_Index SWEAT_Index
#> 3.750147453 44.354838710 106.416793661
#> STP_fix STP_new STP_fix_LM
#> 0.359981513 0.200492750 0.127176069
#> STP_new_LM SCP_fix SCP_new
#> 0.070831081 6.233760868 4.924259548
#> SCP_fix_LM SCP_new_LM SHIP
#> -0.536673003 -0.423936243 0.628661858
#> HSI DCP MU_WMAXSHEAR
#> 1.715903102 1.150670745 1213.184834430
#> SB_WMAXSHEAR ML_WMAXSHEAR MU_EFF_WMAXSHEAR
#> 1213.184834430 1033.105080761 958.335927637
#> SB_EFF_WMAXSHEAR ML_EFF_WMAXSHEAR EHI_500m
#> 958.335927637 797.354807442 0.279777567
#> EHI_01km EHI_03km EHI_500m_LM
#> 0.562298052 2.163834321 0.098841218
#> EHI_01km_LM EHI_03km_LM SHERBS3
#> 0.198651467 -0.186287458 0.648168875
#> SHERBE SHERBS3_v2 SHERBE_v2
#> 0.701454223 0.864225167 0.935272297
#> DEI DEI_eff TIP
#> 1.519761039 1.188457460 2.435630415Example 4
Download sounding and draw hodograph:
data("northplatte")
sounding_hodograph(ws = northplatte$ws, wd = northplatte$wd,
altitude = northplatte$altitude,max_speed = 40)
title("North Platte - 03.07.1999, 00:00 UTC")
Example 5 - thunder in Python
It is possible to launch thunder under Python via
rpy2 library. Below you can find the minimum reproducible
example. Make sure that pandas and rpy2
libraries are available for your Python environment. If not install
required python packages with pip install pandas and
pip install rpy2
Launch thunder under Python with rpy2 and
numpy:
### load required packages
#> from rpy2.robjects.packages import importr
#> from rpy2.robjects import r,pandas2ri
#> import rpy2.robjects as robjects
#> pandas2ri.activate()
### load thunder package (make sure that it was installed in R before)
#> importr('thunder')
### download North Platte sounding
#> profile = robjects.r['get_sounding'](wmo_id = 72562, yy = 1999, mm = 7, dd = 3,hh = 0)
### compute convective parameters
#> parameters = robjects.r['sounding_compute'](profile['pressure'], profile['altitude'], #> profile['temp'], profile['dpt'], profile['wd'], profile['ws'], accuracy = 2)
### customize output and print all computed variables, e.g. most-unstable CAPE (first element) equals 9413 J/kg
#> print(list(map('{:.2f}'.format, parameters)))
# '9413.29', '233.35', '1713.74', '0.00', '775.00', '775.00', '15500.00', '-16.55',
# '137.21', '-66.63', '23.98', '23.98', '23.36', '9413.29', '233.35', '1713.74', '0.00',
# '775.00', '775.00', '15500.00', '-16.55', '137.21', '-66.63', '23.98', '23.98', '23.36',
# '7805.13', '115.22', '1515.81', '-4.35', '950.00', '950.00', '15000.00', '-14.66',
# '124.94', '-68.41', '22.46', '22.46', '21.17', '-9.57', '-6.68', '-8.80', '-8.68',
# '-9.06', '-7.70', '4250.00', '3500.00', '0.00', '2866.00', '50.57', '52.93', '1381.81',
# '308.98', '29.00', '37.59', '87.03', '0.58', '0.40', '0.47', '8.85', '11.21', '13.88',
# '20.28', '29.33', '6.84', '21.70', '28.32', '28.32', '27.17', '17.06', '12.53', '12.53',
# '11.74', '7.09', '6.08', '7.77', '7.69', '19.89', '62.07', '110.06', '156.48', '6.25',
# '7.77', '4.26', '-42.78', '284.67', '5.65', '197.60', '14.19', '218.89', '7.77',
# '31.50', '-12.14', '60.40', '677.12', '4.67', '6.10', '29.46', '29.46', '3.86', '12.35',
# '2783.07', '2783.07', '2534.22', '3886.07', '3886.07', '3395.00']Important notes
- Remember to always input wind speed data in knots.
- Script will always consider first height level as the surface (h = 0), therefore input height data can be as above sea level (ASL) or above ground level (AGL).
- For efficiency purposes it is highly recommended to input data for a maximum of 16-18 km AGL or lower.
- Values of parameters will be different for different accuracy settings.