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
try:
import seaborn as sns
sns.set(rc={"figure.figsize": (12, 6)})
except ImportError:
print('We suggest you install seaborn using conda or pip and rerun this cell')
# built in python modules
import datetime
import os
# python add-ons
import numpy as np
import pandas as pd
try:
import netCDF4
from netCDF4 import num2date
except ImportError:
print('We suggest you install netCDF4 using conda rerun this cell')
# 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 = 'US/Arizona'
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 [4]:
from pvlib.forecast import GFS, HRRR_ESRL, NAM, NDFD, HRRR, RAP
In [5]:
# GFS model, defaults to 0.5 degree resolution
fm = GFS()
In [6]:
# retrieve data
data = fm.get_data(latitude, longitude, start, end)
In [7]:
data
Out[7]:
In [8]:
data = fm.process_data(data)
In [9]:
data[['ghi', 'dni', 'dhi']].plot()
Out[9]:
In [10]:
cs = fm.location.get_clearsky(data.index)
In [11]:
fig, ax = plt.subplots()
cs['ghi'].plot(ax=ax, label='ineichen')
data['ghi'].plot(ax=ax, label='gfs+liujordan')
ax.set_ylabel('ghi')
ax.legend()
Out[11]:
In [12]:
fig, ax = plt.subplots()
cs['dni'].plot(ax=ax, label='ineichen')
data['dni'].plot(ax=ax, label='gfs+liujordan')
ax.set_ylabel('ghi')
ax.legend()
Out[12]:
In [13]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [14]:
data
Out[14]:
In [15]:
data['temp_air'].plot()
plt.ylabel('temperature (%s)' % fm.units['temp_air'])
Out[15]:
In [16]:
cloud_vars = ['total_clouds', 'low_clouds', 'mid_clouds', 'high_clouds']
In [17]:
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))
Out[17]:
In [18]:
total_cloud_cover = data['total_clouds']
In [19]:
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')
Out[19]:
In [20]:
# GFS model at 0.25 degree resolution
fm = GFS(resolution='quarter')
In [22]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [23]:
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))
Out[23]:
In [24]:
data
Out[24]:
In [25]:
fm = NAM()
In [26]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [27]:
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))
Out[27]:
In [28]:
data['ghi'].plot(linewidth=2, ls='-')
plt.ylabel('GHI W/m**2')
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
Out[28]:
In [29]:
data
Out[29]:
In [30]:
fm = NDFD()
In [31]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [33]:
total_cloud_cover = data['total_clouds']
temp = data['temp_air']
wind = data['wind_speed']
In [34]:
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)
Out[34]:
In [36]:
temp.plot(color='r', linewidth=2)
plt.ylabel('Temperature' + ' (%s)' % fm.units['temp_air'])
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
Out[36]:
In [37]:
wind.plot(color='r', linewidth=2)
plt.ylabel('Wind Speed' + ' (%s)' % fm.units['wind_speed'])
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
Out[37]:
In [38]:
data
Out[38]:
In [39]:
fm = RAP(resolution=20)
In [40]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [41]:
cloud_vars = ['total_clouds', 'high_clouds', 'mid_clouds', 'low_clouds']
In [42]:
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))
Out[42]:
In [43]:
data
Out[43]:
In [44]:
fm = HRRR()
In [46]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [47]:
cloud_vars = ['total_clouds', 'high_clouds', 'mid_clouds', 'low_clouds']
In [48]:
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')
plt.legend(bbox_to_anchor=(1.18,1.0))
Out[48]:
In [49]:
data
Out[49]:
In [50]:
fm = HRRR_ESRL()
In [51]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [ ]:
cloud_vars = ['total_clouds','high_clouds','mid_clouds','low_clouds']
In [ ]:
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 [ ]:
data['ghi'].plot(linewidth=2, ls='-')
plt.ylabel('GHI W/m**2')
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
In [54]:
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__CEC_2014_']
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
irradiance = fx_data[['ghi', 'dni', 'dhi']]
weather = fx_data[['wind_speed', 'temp_air']]
mc.run_model(fx_data.index, irradiance=irradiance, weather=weather)
Out[54]:
In [55]:
mc.total_irrad.plot()
Out[55]:
In [56]:
mc.temps.plot()
Out[56]:
In [57]:
mc.ac.plot()
Out[57]:
In [ ]: