aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorGravatar Jonas Gunz <himself@jonasgunz.de> 2023-08-24 17:22:32 +0200
committerGravatar Jonas Gunz <himself@jonasgunz.de> 2023-08-24 17:22:32 +0200
commitb85a2929e76ddc39f83ac0403f8356e05b71d129 (patch)
tree17a8e630d55d1d95bd6f071b07851665349e8586
parent60ab2d134df00f8099add774e698b2d7d2395bdc (diff)
downloadmeteo_toolbox-b85a2929e76ddc39f83ac0403f8356e05b71d129.tar.gz
stuff
-rw-r--r--Readme.md5
-rwxr-xr-xgrib.py23
-rwxr-xr-xicon_download.sh17
-rw-r--r--requirements.txt2
-rw-r--r--skewt.py70
-rwxr-xr-xxar.py34
6 files changed, 141 insertions, 10 deletions
diff --git a/Readme.md b/Readme.md
index 32de887..38ed16c 100644
--- a/Readme.md
+++ b/Readme.md
@@ -3,3 +3,8 @@
DTKE_CON for convection forecast?
w_ctmax updraft velocity
+
+## Sources
+
+DWD OpenData: opendata.dwd.de
+ECMWF OpenData: https://www.ecmwf.int/en/forecasts/datasets/open-data
diff --git a/grib.py b/grib.py
index c298450..767f3d7 100755
--- a/grib.py
+++ b/grib.py
@@ -16,15 +16,34 @@ target_lon='10.0646952'
grib = pygrib.open('dwd_icon-d2/combined.grib2')
+# GRIB-Objekt wichtige attribute:
+# analDate
+# validDate
+# units
+# level
+# typeOfLevel
+# name
+# shortName
for grb in grib:
- print(grb)
+ #print(grb)
vals = grb.values
lats, lons = grb.latlons()
+ for k in grb.keys():
+ print(k)
+
+ print(grb.analDate)
+ print(grb.validDate)
+ print(grb.units)
+ print(grb.typeOfLevel)
+ break
+
+ # TODO Data in xarray?
+
print(grb.name)
- #print(grb.level)
+ print(grb.level)
#print('lats min/max: ', lats.shape, lats.max(), lats.min())
#print('lons min/max: ', lons.shape, lons.max(), lons.min())
diff --git a/icon_download.sh b/icon_download.sh
index cf00cca..343b10c 100755
--- a/icon_download.sh
+++ b/icon_download.sh
@@ -9,9 +9,10 @@ MODEL=icon-d2
MODEL_LONG=icon-d2_germany
BASE="http://opendata.dwd.de/weather/nwp"
-RUN="15"
+RUN="00"
PARAMETERS=( "t" "relhum" "u" "v" "fi" )
-PARAMETERS_SINGLE_LEVEL=( "cape_ml" "cin_ml" "tot_prec" "w_ctmax" )
+# tot_prec and cape_ml/cin_ml is in 15min intervals and screws with xygrib
+PARAMETERS_SINGLE_LEVEL=( "w_ctmax" )
PRESSURE_LEVELS=( "1000" "975" "950" "850" "700" "600" "500" "400" "300" "250" "200" )
OFFSETS=( "000" "003" "006" "009" "012" "015" "018" "024" )
DATE=$(date +%Y%m%d)
@@ -20,24 +21,24 @@ mkdir -p $OUTDIR
echo -n > "$OUTDIR/index.txt"
-for OFFSET in ${OFFSETS[@]}; do
- for PARAMETER in ${PARAMETERS[@]}; do
- for LEVEL in ${PRESSURE_LEVELS[@]}; do
+for OFFSET in "${OFFSETS[@]}"; do
+ for PARAMETER in "${PARAMETERS[@]}"; do
+ for LEVEL in "${PRESSURE_LEVELS[@]}"; do
URL="$BASE/$MODEL/grib/$RUN/$PARAMETER/${MODEL_LONG}_regular-lat-lon_pressure-level_${DATE}${RUN}_${OFFSET}_${LEVEL}_${PARAMETER}.grib2.bz2"
BNAME=$(basename "$URL")
echo Getting "$URL"
echo "${BNAME%.bz2}" >> $OUTDIR/index.txt
- wget -q --directory-prefix=$OUTDIR "$URL"
+ wget -q --directory-prefix=$OUTDIR "$URL" || echo FAILED!
done
done
- for PARAMETER in ${PARAMETERS_SINGLE_LEVEL[@]}; do
+ for PARAMETER in "${PARAMETERS_SINGLE_LEVEL[@]}"; do
URL="$BASE/$MODEL/grib/$RUN/$PARAMETER/${MODEL_LONG}_regular-lat-lon_single-level_${DATE}${RUN}_${OFFSET}_2d_${PARAMETER}.grib2.bz2"
BNAME=$(basename "$URL")
echo Getting "$URL"
echo "${BNAME%.bz2}" >> $OUTDIR/index.txt
- wget -q --directory-prefix=$OUTDIR "$URL"
+ wget -q --directory-prefix=$OUTDIR "$URL" || echo FAILED!
done
done
diff --git a/requirements.txt b/requirements.txt
index 7ce482f..aef6541 100644
--- a/requirements.txt
+++ b/requirements.txt
@@ -1,2 +1,4 @@
metpy
pygrib
+xarray
+cfgrib
diff --git a/skewt.py b/skewt.py
new file mode 100644
index 0000000..43bcc91
--- /dev/null
+++ b/skewt.py
@@ -0,0 +1,70 @@
+import matplotlib.gridspec as gridspec
+import matplotlib.pyplot as plt
+
+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
+
+class Skewt:
+ def __init__(self, p, T, Td, title=None):
+ self._p = p
+ self._T = T
+ self._Td = Td
+
+
+ # Create a new figure. The dimensions here give a good aspect ratio
+ self._fig = plt.figure(figsize=(9, 9))
+ plt.rcParams["font.family"] = "monospace"
+ #self._fig = plt.figure()
+
+ if title is not None:
+ plt.suptitle(title, x=0, y=0, va='bottom', ha='left')
+
+ # Grid for plots
+ self._gs = gridspec.GridSpec(3, 3)
+ self._skew = SkewT(self._fig, rotation=45, subplot=self._gs[:, :2])
+
+ # Plot the data using normal plotting functions, in this case using
+ # log scaling in Y, as dictated by the typical meteorological plot
+ self._skew.plot(p, T, 'r')
+ self._skew.plot(p, Td, 'b')
+
+ plt.xlabel('$T$ $[^o C]$')
+ plt.ylabel('$p$ $[hPa]$')
+
+ def addWindUV(self, u, v):
+ self._u = u
+ self._v = v
+
+ ax = self._fig.add_subplot(self._gs[0, -1])
+ h = Hodograph(ax, component_range=max(u + v).magnitude)
+ h.add_grid(increment=20)
+ h.plot_colormapped(u, v, self._p)
+ plt.tight_layout()
+ plt.xlabel('$m/s$')
+ plt.ylabel('$m/s$')
+ self._skew.plot_barbs(self._p, u, v)
+
+ def addInfo(self, lines):
+ # TODO
+ ax = self._fig.add_subplot(self._gs[1,-1])
+ ax.text(0, 0, lines, ha='left', va='center', size=12)
+ ax.axis("off")
+
+
+ def plot(self, filename=None):
+ # Add the relevant special lines
+ #self._skew.ax.set_ylim(max(self._p), min(self._p))
+ self._skew.ax.set_ylim(1000, 100)
+ self._skew.plot_dry_adiabats()
+ self._skew.plot_moist_adiabats()
+ self._skew.plot_mixing_lines()
+
+ # Good bounds for aspect ratio
+ self._skew.ax.set_xlim(-30, 40)
+
+ if filename is not None:
+ plt.savefig(filename)
+ else:
+ plt.show()
diff --git a/xar.py b/xar.py
new file mode 100755
index 0000000..fca3887
--- /dev/null
+++ b/xar.py
@@ -0,0 +1,34 @@
+#!/usr/bin/env python3
+
+import xarray as xr
+from metpy.units import units
+import metpy.calc as mpcalc
+
+import skewt
+
+#grib = pygrib.open('dwd_icon-d2/combined.grib2')
+data = xr.load_dataset('dwd_icon-d2/combined.grib2', engine='cfgrib')
+
+lat = 47.9626
+lon = 11.9964
+
+for_temp = data.sel(latitude=lat, longitude = lon, method='nearest')
+
+for_temp = for_temp[['r', 't', 'u', 'v']]
+
+for step in for_temp.coords['step']:
+ this_step = for_temp.sel(step=step)
+
+ p = this_step.coords['isobaricInhPa'].values * units.hPa
+ T = this_step.t.values * units.K
+ relHum = this_step.r.values * units.percent
+ Td = mpcalc.dewpoint_from_relative_humidity(T, relHum)
+ u = this_step.u.values * (units.m / units.s)
+ v = this_step.v.values * (units.m / units.s)
+
+ skt = skewt.Skewt(p=p, T=T, Td=Td)
+ skt.addWindUV(u, v)
+ skt.addInfo("TEST")
+ skt.plot(filename=f'skewt.png')
+
+ break