1# Copyright (c) 2017 MetPy Developers.
2# Distributed under the terms of the BSD 3-Clause License.
3# SPDX-License-Identifier: BSD-3-Clause
4"""
5Meteogram
6=========
7
8Plots time series data as a meteogram.
9"""
10
11import datetime as dt
12
13import matplotlib as mpl
14import matplotlib.pyplot as plt
15import numpy as np
16
17from metpy.calc import dewpoint_from_relative_humidity
18from metpy.cbook import get_test_data
19from metpy.plots import add_metpy_logo
20from metpy.units import units
21
22
23def calc_mslp(t, p, h):
24    return p * (1 - (0.0065 * h) / (t + 0.0065 * h + 273.15)) ** (-5.257)
25
26
27# Make meteogram plot
28class Meteogram:
29    """ Plot a time series of meteorological data from a particular station as a
30    meteogram with standard variables to visualize, including thermodynamic,
31    kinematic, and pressure. The functions below control the plotting of each
32    variable.
33    TO DO: Make the subplot creation dynamic so the number of rows is not
34    static as it is currently. """
35
36    def __init__(self, fig, dates, probeid, time=None, axis=0):
37        """
38        Required input:
39            fig: figure object
40            dates: array of dates corresponding to the data
41            probeid: ID of the station
42        Optional Input:
43            time: Time the data is to be plotted
44            axis: number that controls the new axis to be plotted (FOR FUTURE)
45        """
46        if not time:
47            time = dt.datetime.utcnow()
48        self.start = dates[0]
49        self.fig = fig
50        self.end = dates[-1]
51        self.axis_num = 0
52        self.dates = mpl.dates.date2num(dates)
53        self.time = time.strftime('%Y-%m-%d %H:%M UTC')
54        self.title = f'Latest Ob Time: {self.time}\nProbe ID: {probeid}'
55
56    def plot_winds(self, ws, wd, wsmax, plot_range=None):
57        """
58        Required input:
59            ws: Wind speeds (knots)
60            wd: Wind direction (degrees)
61            wsmax: Wind gust (knots)
62        Optional Input:
63            plot_range: Data range for making figure (list of (min,max,step))
64        """
65        # PLOT WIND SPEED AND WIND DIRECTION
66        self.ax1 = fig.add_subplot(4, 1, 1)
67        ln1 = self.ax1.plot(self.dates, ws, label='Wind Speed')
68        self.ax1.fill_between(self.dates, ws, 0)
69        self.ax1.set_xlim(self.start, self.end)
70        if not plot_range:
71            plot_range = [0, 20, 1]
72        self.ax1.set_ylabel('Wind Speed (knots)', multialignment='center')
73        self.ax1.set_ylim(plot_range[0], plot_range[1], plot_range[2])
74        self.ax1.grid(b=True, which='major', axis='y', color='k', linestyle='--',
75                      linewidth=0.5)
76        ln2 = self.ax1.plot(self.dates, wsmax, '.r', label='3-sec Wind Speed Max')
77
78        ax7 = self.ax1.twinx()
79        ln3 = ax7.plot(self.dates, wd, '.k', linewidth=0.5, label='Wind Direction')
80        ax7.set_ylabel('Wind\nDirection\n(degrees)', multialignment='center')
81        ax7.set_ylim(0, 360)
82        ax7.set_yticks(np.arange(45, 405, 90))
83        ax7.set_yticklabels(['NE', 'SE', 'SW', 'NW'])
84        lines = ln1 + ln2 + ln3
85        labs = [line.get_label() for line in lines]
86        ax7.xaxis.set_major_formatter(mpl.dates.DateFormatter('%d/%H UTC'))
87        ax7.legend(lines, labs, loc='upper center',
88                   bbox_to_anchor=(0.5, 1.2), ncol=3, prop={'size': 12})
89
90    def plot_thermo(self, t, td, plot_range=None):
91        """
92        Required input:
93            T: Temperature (deg F)
94            TD: Dewpoint (deg F)
95        Optional Input:
96            plot_range: Data range for making figure (list of (min,max,step))
97        """
98        # PLOT TEMPERATURE AND DEWPOINT
99        if not plot_range:
100            plot_range = [10, 90, 2]
101        self.ax2 = fig.add_subplot(4, 1, 2, sharex=self.ax1)
102        ln4 = self.ax2.plot(self.dates, t, 'r-', label='Temperature')
103        self.ax2.fill_between(self.dates, t, td, color='r')
104
105        self.ax2.set_ylabel('Temperature\n(F)', multialignment='center')
106        self.ax2.grid(b=True, which='major', axis='y', color='k', linestyle='--',
107                      linewidth=0.5)
108        self.ax2.set_ylim(plot_range[0], plot_range[1], plot_range[2])
109
110        ln5 = self.ax2.plot(self.dates, td, 'g-', label='Dewpoint')
111        self.ax2.fill_between(self.dates, td, self.ax2.get_ylim()[0], color='g')
112
113        ax_twin = self.ax2.twinx()
114        ax_twin.set_ylim(plot_range[0], plot_range[1], plot_range[2])
115        lines = ln4 + ln5
116        labs = [line.get_label() for line in lines]
117        ax_twin.xaxis.set_major_formatter(mpl.dates.DateFormatter('%d/%H UTC'))
118
119        self.ax2.legend(lines, labs, loc='upper center',
120                        bbox_to_anchor=(0.5, 1.2), ncol=2, prop={'size': 12})
121
122    def plot_rh(self, rh, plot_range=None):
123        """
124        Required input:
125            RH: Relative humidity (%)
126        Optional Input:
127            plot_range: Data range for making figure (list of (min,max,step))
128        """
129        # PLOT RELATIVE HUMIDITY
130        if not plot_range:
131            plot_range = [0, 100, 4]
132        self.ax3 = fig.add_subplot(4, 1, 3, sharex=self.ax1)
133        self.ax3.plot(self.dates, rh, 'g-', label='Relative Humidity')
134        self.ax3.legend(loc='upper center', bbox_to_anchor=(0.5, 1.22), prop={'size': 12})
135        self.ax3.grid(b=True, which='major', axis='y', color='k', linestyle='--',
136                      linewidth=0.5)
137        self.ax3.set_ylim(plot_range[0], plot_range[1], plot_range[2])
138
139        self.ax3.fill_between(self.dates, rh, self.ax3.get_ylim()[0], color='g')
140        self.ax3.set_ylabel('Relative Humidity\n(%)', multialignment='center')
141        self.ax3.xaxis.set_major_formatter(mpl.dates.DateFormatter('%d/%H UTC'))
142        axtwin = self.ax3.twinx()
143        axtwin.set_ylim(plot_range[0], plot_range[1], plot_range[2])
144
145    def plot_pressure(self, p, plot_range=None):
146        """
147        Required input:
148            P: Mean Sea Level Pressure (hPa)
149        Optional Input:
150            plot_range: Data range for making figure (list of (min,max,step))
151        """
152        # PLOT PRESSURE
153        if not plot_range:
154            plot_range = [970, 1030, 2]
155        self.ax4 = fig.add_subplot(4, 1, 4, sharex=self.ax1)
156        self.ax4.plot(self.dates, p, 'm', label='Mean Sea Level Pressure')
157        self.ax4.set_ylabel('Mean Sea\nLevel Pressure\n(mb)', multialignment='center')
158        self.ax4.set_ylim(plot_range[0], plot_range[1], plot_range[2])
159
160        axtwin = self.ax4.twinx()
161        axtwin.set_ylim(plot_range[0], plot_range[1], plot_range[2])
162        axtwin.fill_between(self.dates, p, axtwin.get_ylim()[0], color='m')
163        axtwin.xaxis.set_major_formatter(mpl.dates.DateFormatter('%d/%H UTC'))
164
165        self.ax4.legend(loc='upper center', bbox_to_anchor=(0.5, 1.2), prop={'size': 12})
166        self.ax4.grid(b=True, which='major', axis='y', color='k', linestyle='--',
167                      linewidth=0.5)
168        # OTHER OPTIONAL AXES TO PLOT
169        # plot_irradiance
170        # plot_precipitation
171
172
173# set the starttime and endtime for plotting, 24 hour range
174endtime = dt.datetime(2016, 3, 31, 22, 0, 0, 0)
175starttime = endtime - dt.timedelta(hours=24)
176
177# Height of the station to calculate MSLP
178hgt_example = 292.
179
180
181# Parse dates from .csv file, knowing their format as a string and convert to datetime
182def parse_date(date):
183    return dt.datetime.strptime(date.decode('ascii'), '%Y-%m-%d %H:%M:%S')
184
185
186testdata = np.genfromtxt(get_test_data('timeseries.csv', False), names=True, dtype=None,
187                         usecols=list(range(1, 8)),
188                         converters={'DATE': parse_date}, delimiter=',')
189
190# Temporary variables for ease
191temp = testdata['T']
192pressure = testdata['P']
193rh = testdata['RH']
194ws = testdata['WS']
195wsmax = testdata['WSMAX']
196wd = testdata['WD']
197date = testdata['DATE']
198
199# ID For Plotting on Meteogram
200probe_id = '0102A'
201
202data = {'wind_speed': (np.array(ws) * units('m/s')).to(units('knots')),
203        'wind_speed_max': (np.array(wsmax) * units('m/s')).to(units('knots')),
204        'wind_direction': np.array(wd) * units('degrees'),
205        'dewpoint': dewpoint_from_relative_humidity((np.array(temp) * units.degC).to(units.K),
206                                                    np.array(rh) / 100.).to(units('degF')),
207        'air_temperature': (np.array(temp) * units('degC')).to(units('degF')),
208        'mean_slp': calc_mslp(np.array(temp), np.array(pressure), hgt_example) * units('hPa'),
209        'relative_humidity': np.array(rh), 'times': np.array(date)}
210
211fig = plt.figure(figsize=(20, 16))
212add_metpy_logo(fig, 250, 180)
213meteogram = Meteogram(fig, data['times'], probe_id)
214meteogram.plot_winds(data['wind_speed'], data['wind_direction'], data['wind_speed_max'])
215meteogram.plot_thermo(data['air_temperature'], data['dewpoint'])
216meteogram.plot_rh(data['relative_humidity'])
217meteogram.plot_pressure(data['mean_slp'])
218fig.subplots_adjust(hspace=0.5)
219plt.show()
220