Forecast Tutorial¶
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:
- Python 3.5.2
- IPython 5.0.0
- pandas 0.18.0
- matplotlib 1.5.1
- netcdf4 1.2.1
- siphon 0.4.0
It should work with other Python and Pandas versions. It requires pvlib >= 0.3.0 and IPython >= 3.0.
Authors:
- Derek Groenendyk (@moonraker), University of Arizona, November 2015
- Will Holmgren (@wholmgren), University of Arizona, November 2015, January 2016, April 2016, July 2016
Setup¶
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
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)
GFS (0.5 deg)¶
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
Out[6]:
In [7]:
data = fm.process_data(data)
In [8]:
data[['ghi', 'dni', 'dhi']].plot()
Out[8]:
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+liujordan')
ax.set_ylabel('ghi')
ax.legend()
Out[10]:
In [11]:
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[11]:
In [12]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [13]:
data
Out[13]:
In [14]:
data['temp_air'].plot()
plt.ylabel('temperature (%s)' % fm.units['temp_air'])
Out[14]:
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))
Out[16]:
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')
Out[18]:
GFS (0.25 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))
Out[21]:
In [22]:
data
Out[22]:
NAM¶
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))
Out[25]:
In [26]:
data['ghi'].plot(linewidth=2, ls='-')
plt.ylabel('GHI W/m**2')
plt.xlabel('Forecast Time ('+str(data.index.tz)+')')
Out[26]:
In [27]:
data
Out[27]:
NDFD¶
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)
Out[31]:
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
Out[34]:
RAP¶
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))
Out[38]:
In [39]:
data
Out[39]:
HRRR¶
In [40]:
fm = HRRR()
In [41]:
# retrieve data
data = fm.get_processed_data(latitude, longitude, start, end)
In [42]:
cloud_vars = ['total_clouds', 'high_clouds', 'mid_clouds', 'low_clouds']
In [43]:
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[43]:
In [44]:
data
Out[44]:
HRRR (ESRL)¶
In [45]:
fm = HRRR_ESRL()
In [46]:
# 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)+')')
Quick power calculation¶
In [48]:
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
mc.run_model(fx_data.index, weather=fx_data)
Out[48]:
In [49]:
mc.total_irrad.plot()
Out[49]:
In [50]:
mc.temps.plot()
Out[50]:
In [51]:
mc.ac.plot()
Out[51]:
In [ ]: