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
106
107
|
#!/usr/bin/env python3
# https://unidata.github.io/MetPy/latest/examples/plots/Skew-T_Layout.html#sphx-glr-examples-plots-skew-t-layout-py
import matplotlib.gridspec as gridspec
import matplotlib.pyplot as plt
import csv
import metpy.calc as mpcalc
from metpy.cbook import get_test_data
from metpy.plots import add_metpy_logo, Hodograph, SkewT
from metpy.units import units
import datetime
import requests
STATION=10548
DATE=(datetime.date.today() - datetime.timedelta(days = 1)).strftime('%Y-%m-%d')
TIME='23:00:00'
# http://weather.uwyo.edu/cgi-bin/bufrraob.py?datetime=2023-07-20%2023:00:00&id=10548&type=TEXT:CSV
URL=f'http://weather.uwyo.edu/cgi-bin/bufrraob.py?datetime={DATE}%20{TIME}&id={STATION}&type=TEXT:CSV'
print(URL)
result = requests.get(URL)
RAW_DATA_LINES=result.content.decode('UTF-8').splitlines()
#RAW_DATA_LINES= open('data.csv').readlines()
p = []
T = []
Td = []
wind_speed = []
wind_dir = []
csvreader = csv.reader(RAW_DATA_LINES, dialect='excel')
next(csvreader) # Skip header
for row in csvreader:
if sum(map(lambda s : len(s.strip()) == 0, row)):
# skip any line with empty values
continue
if float(row[3]) in p: # Skip double p entries
continue
p.append(float(row[3]))
T.append(float(row[5]))
Td.append(float(row[6]))
wind_speed.append(float(row[12]))
wind_dir.append(float(row[11]))
p = p * units.hPa
T = T * units.degC
Td = Td * units.degC
wind_speed = wind_speed * units.knots
wind_dir = wind_dir * units.degrees
u, v = mpcalc.wind_components(wind_speed, wind_dir)
ground_parcel= mpcalc.parcel_profile(p, T[0], Td[0])
ccl = mpcalc.ccl(p,T,Td)
lcl = mpcalc.lcl(p[0],T[0],Td[0])
ccl_ground=mpcalc.dry_lapse(p[:1], ccl[1], ccl[0])
ccl_ground_parcel= mpcalc.parcel_profile(p, ccl_ground[0], Td[0])
# Create a new figure. The dimensions here give a good aspect ratio
fig = plt.figure(figsize=(9, 9))
# Grid for plots
gs = gridspec.GridSpec(3, 3)
skew = SkewT(fig, rotation=45, subplot=gs[:, :])
# Plot the data using normal plotting functions, in this case using
# log scaling in Y, as dictated by the typical meteorological plot
skew.plot(p, T, 'r')
skew.plot(p, Td, 'g')
skew.plot(p, ground_parcel, 'b')
skew.plot(p, ccl_ground_parcel, 'y')
barb_div = int(len(p)/20)
skew.plot_barbs(p[::barb_div], u[::barb_div], v[::barb_div])
skew.plot(lcl[0], lcl[1], 'o', markerfacecolor='black')
skew.plot(ccl[0], ccl[1], 'o', markerfacecolor='red')
skew.shade_cape(p,T,ground_parcel)
skew.shade_cin(p,T,ground_parcel)
skew.ax.set_ylim(1000, 100)
# Add the relevant special lines
skew.plot_dry_adiabats()
skew.plot_moist_adiabats()
skew.plot_mixing_lines()
# Good bounds for aspect ratio
skew.ax.set_xlim(-30, 40)
# Create a hodograph
ax = fig.add_subplot(gs[0, -1])
h = Hodograph(ax, component_range=60.)
h.add_grid(increment=20)
h.plot(u[::barb_div], v[::barb_div])
# Show the plot
plt.show()
|