This tutorial will walk through forecast data from Unidata forecast model data using the forecast.py module within pvlib.
Table of contents:
This tutorial has been tested against the following package versions:
It should work with other Python and Pandas versions. It requires pvlib >= 0.3.0 and IPython >= 3.0.
Authors:
In [1]:
%matplotlib inline
import matplotlib.pyplot as plt
# built in python modules
import datetime
import os
# python add-ons
import numpy as np
import pandas as pd
# for accessing UNIDATA THREDD servers
from siphon.catalog import TDSCatalog
from siphon.ncss import NCSS
import pvlib
from pvlib.forecast import GFS, HRRR_ESRL, NAM, NDFD, HRRR, RAP
In [2]:
# Choose a location and time.
# Tucson, AZ
latitude = 32.2
longitude = -110.9
tz = 'America/Phoenix'
start = pd.Timestamp(datetime.date.today(), tz=tz) # today's date
end = start + pd.Timedelta(days=7) # 7 days from today
print(start, end)
In [3]:
from pvlib.forecast import GFS, HRRR_ESRL, NAM, NDFD, HRRR, RAP
In [4]:
# GFS model, defaults to 0.5 degree resolution
fm = GFS()
In [5]:
# retrieve data
data = fm.get_data(latitude, longitude, start, end)
In [6]:
data[sorted(data.columns)]
Out[6]:
In [7]:
data = fm.process_data(data)
In [8]:
data[['ghi', 'dni', 'dhi']].plot();
In [9]:
cs = fm.location.get_clearsky(data.index)
In [10]:
fig, ax = plt.subplots()
cs['ghi'].plot(ax=ax, label='ineichen')
data['ghi'].plot(ax=ax, label='gfs+larson')
ax.set_ylabel('ghi')
ax.legend();
In [11]:
fig, ax = plt.subplots()
cs['dni'].plot(ax=ax, label='ineichen')
data['dni'].plot(ax=ax, label='gfs+larson')
ax.set_ylabel('ghi')
ax.legend();
In [12]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [13]:
data[sorted(data.columns)]
Out[13]:
In [14]:
data['temp_air'].plot()
plt.ylabel('temperature (%s)' % fm.units['temp_air']);
In [15]:
cloud_vars = ['total_clouds', 'low_clouds', 'mid_clouds', 'high_clouds']
In [16]:
for varname in cloud_vars:
data[varname].plot()
plt.ylabel('Cloud cover' + ' %')
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
plt.title('GFS 0.5 deg')
plt.legend(bbox_to_anchor=(1.18,1.0));
In [17]:
total_cloud_cover = data['total_clouds']
In [18]:
total_cloud_cover.plot(color='r', linewidth=2)
plt.ylabel('Total cloud cover' + ' (%s)' % fm.units['total_clouds'])
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
plt.title('GFS 0.5 deg');
In [19]:
# GFS model at 0.25 degree resolution
fm = GFS(resolution='quarter')
In [20]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [21]:
for varname in cloud_vars:
data[varname].plot(ls='-', linewidth=2)
plt.ylabel('Cloud cover' + ' %')
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
plt.title('GFS 0.25 deg')
plt.legend(bbox_to_anchor=(1.18,1.0));
In [22]:
data[sorted(data.columns)]
Out[22]:
In [23]:
fm = NAM()
In [24]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [25]:
for varname in cloud_vars:
data[varname].plot(ls='-', linewidth=2)
plt.ylabel('Cloud cover' + ' %')
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
plt.title('NAM')
plt.legend(bbox_to_anchor=(1.18,1.0));
In [26]:
data['ghi'].plot(linewidth=2, ls='-')
plt.ylabel('GHI W/m**2')
plt.xlabel('Forecast Time ('+str(data.index.tz)+')');
In [27]:
data[sorted(data.columns)]
Out[27]:
In [28]:
fm = NDFD()
In [29]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [30]:
total_cloud_cover = data['total_clouds']
temp = data['temp_air']
wind = data['wind_speed']
In [31]:
total_cloud_cover.plot(color='r', linewidth=2)
plt.ylabel('Total cloud cover' + ' (%s)' % fm.units['total_clouds'])
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
plt.title('NDFD')
plt.ylim(0,100);
In [32]:
temp.plot(color='r', linewidth=2)
plt.ylabel('Temperature' + ' (%s)' % fm.units['temp_air'])
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
Out[32]:
In [33]:
wind.plot(color='r', linewidth=2)
plt.ylabel('Wind Speed' + ' (%s)' % fm.units['wind_speed'])
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
Out[33]:
In [34]:
data[sorted(data.columns)]
Out[34]:
In [35]:
fm = RAP(resolution=20)
In [36]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [37]:
cloud_vars = ['total_clouds', 'high_clouds', 'mid_clouds', 'low_clouds']
In [38]:
for varname in cloud_vars:
data[varname].plot(ls='-', linewidth=2)
plt.ylabel('Cloud cover' + ' %')
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
plt.title('RAP')
plt.legend(bbox_to_anchor=(1.18,1.0));
In [39]:
data[sorted(data.columns)]
Out[39]:
In [40]:
fm = HRRR()
In [41]:
data_raw = fm.get_data(latitude, longitude, start, end)
In [42]:
# The HRRR model pulls in u, v winds for 2 layers above ground (10 m, 80 m)
# They are labeled as _0, _1 in the raw data
data_raw[sorted(data_raw.columns)]
Out[42]:
In [43]:
data = fm.get_processed_data(latitude, longitude, start, end)
In [44]:
cloud_vars = ['total_clouds', 'high_clouds', 'mid_clouds', 'low_clouds']
In [45]:
for varname in cloud_vars:
data[varname].plot(ls='-', linewidth=2)
plt.ylabel('Cloud cover' + ' %')
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
plt.title('RAP')
plt.legend(bbox_to_anchor=(1.18,1.0));
In [46]:
data['temp_air'].plot(color='r', linewidth=2)
plt.ylabel('Temperature' + ' (%s)' % fm.units['temp_air'])
plt.xlabel('Forecast Time ('+str(data.index.tz)+')');
In [47]:
data['wind_speed'].plot(color='r', linewidth=2)
plt.ylabel('Wind Speed' + ' (%s)' % fm.units['wind_speed'])
plt.xlabel('Forecast Time ('+str(data.index.tz)+')');
In [48]:
data[sorted(data.columns)]
Out[48]:
In [49]:
# NBVAL_SKIP
fm = HRRR_ESRL()
In [ ]:
# retrieve data
# NBVAL_SKIP
data = fm.get_processed_data(latitude, longitude, start, end)
In [ ]:
# NBVAL_SKIP
cloud_vars = ['total_clouds','high_clouds','mid_clouds','low_clouds']
In [ ]:
# NBVAL_SKIP
for varname in cloud_vars:
data[varname].plot(ls='-', linewidth=2)
plt.ylabel('Cloud cover' + ' %')
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
plt.title('HRRR_ESRL')
plt.legend(bbox_to_anchor=(1.18,1.0));
In [ ]:
# NBVAL_SKIP
data['ghi'].plot(linewidth=2, ls='-')
plt.ylabel('GHI W/m**2')
plt.xlabel('Forecast Time ('+str(data.index.tz)+')');
In [51]:
from pvlib.pvsystem import PVSystem, retrieve_sam
from pvlib.modelchain import ModelChain
sandia_modules = retrieve_sam('SandiaMod')
sapm_inverters = retrieve_sam('cecinverter')
module = sandia_modules['Canadian_Solar_CS5P_220M___2009_']
inverter = sapm_inverters['ABB__MICRO_0_25_I_OUTD_US_208__208V_']
system = PVSystem(module_parameters=module,
inverter_parameters=inverter)
# fx is a common abbreviation for forecast
fx_model = GFS()
fx_data = fx_model.get_processed_data(latitude, longitude, start, end)
# use a ModelChain object to calculate modeling intermediates
mc = ModelChain(system, fx_model.location,
orientation_strategy='south_at_latitude_tilt')
# extract relevant data for model chain
mc.run_model(weather=fx_data)
Out[51]:
In [52]:
mc.total_irrad.plot();
In [53]:
mc.cell_temperature.plot();
In [54]:
mc.ac.plot();
In [ ]: