Global Horizontal Irradiance data for Canada

Environment Canada publishes Canadian Weather Energy and Engineering Datasets (CWEEDS), which contains a zip file with several years of weather observations for each province. These datasets contain solar radiation information which is very useful for solar panel modeling.

To get the mean hourly values of  Extraterrestrial, Global Horizontal, Direct Normal and Diffuse Horizontal radiation for a given location, here is python code that uses pandas:

import datetime
import os.path
import pandas as pd

def WY2dateParser(str):
    """Parses dates in the WY2 format, which uses hours 1 through 24."""
    try:
        time = datetime.datetime.strptime(str, "%Y%m%d%H")
    except ValueError:
        # replace 24 at the end of str with 23
        str = str[:-2] + '23';
        time = datetime.datetime.strptime(str, "%Y%m%d%H")
        time += datetime.timedelta(hours=1)
    return time


def MonthDayHourTupleToDateTime(MDHtuple, year):
    return datetime.datetime(year, MDHtuple[0], MDHtuple[1], MDHtuple[2], 0)

def hourlyIrradianceData(documentPath):
    """
    :param documentPath: file path of CWEED PY2 document to parse
    :return: pandas dataframe with (month, day, hour) tuple indexes and average
             irradiance values values in kJ/m^2 columns
    """

    if not os.path.isfile(documentPath):
        raise Exception(documentPath + " is not a valid file path.")

    radData = pd.read_fwf(documentPath, colspecs=[(6, 16), (16, 20), (20, 24), (26, 30), (32,36)],
                       names=['timestamp', 'Extraterrestrial', 'Global Horizontal', 'Direct Normal', 'Diffuse Horizontal'],
                       na_values="9999", parse_dates=True, index_col=0, date_parser=WY2dateParser)
    return radData.groupby(lambda x: (x.month, x.day, x.hour)).mean()

To use it, first download the zip file of your province and extract it to your home directory. Then you can use the above code to get for example the hourly Global Horizontal irradiance for June 2015 like this:

from os.path import expanduser
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import axes3d, Axes3D
from matplotlib import cm
%matplotlib inline

numberOfPanels = 20
panelSurface = 1.6434 # in m^2
totalPanelSurface = numberOfPanels * panelSurface

# plot Global Horizontal Irradiance for June 2015
start = datetime.datetime(2015,6,1)
stop = datetime.datetime(2015,6,30)

x,y,z = [],[],[]
radData = hourlyIrradianceData(expanduser('~/ALBERTA/EdmontonNamao_1956-1994/CANA6.WY2'))
# kJ/m^2 to wH
pannelIrradianceData = radData * (1000.0 * totalPanelSurface / 3600.00)

stopNextDay = stop + datetime.timedelta(days=1)
inputData = pannelIrradianceData[(start.month, start.day, 1):(stopNextDay.month, stopNextDay.day, 0)]['Global Horizontal']
inputData.index = inputData.index.map(lambda x: MonthDayHourTupleToDateTime(x, start.year))

for timestamp in inputData.index:
    x.append(timestamp.day)
    y.append(timestamp.hour)
    energy = inputData[timestamp]
    z.append(energy)
        
fig = plt.figure()
fig.set_size_inches(12, 6)
ax = Axes3D(fig) # Use ax = fig.add_subplot(111, projection='3d') for more recent versions of matplotlib
ax.plot_trisurf(x,y,z,cmap=cm.coolwarm, linewidth=0)
ax.set_xlabel('Day')
ax.set_ylabel('Hour')
ax.set_zlabel('Wh')

plt.show()