One part of the Clark et al. (2015) study explored the impact of the choice of canopy shortwave radiation parameterizations on simulations of below canopy shortwave radiation for three representative water years at the aspen site in the Reynolds Mountain East catchment. This study looked at four different canopy shortwave radiation parameterizations: BeersLaw method(as implemented in VIC), NL_scatter method(Nijssen and Lettenmaier, JGR 1999:NL 1999), UEB_2stream method(Mahat and Tarboton, WRR 2011:MT 2012), CLM_2stream method(Dick 1983)

In this Jupyter Notebook, the pySUMMA library is used to reproduce this analysis. First, the four different canopy shortwave radiation parameteriations are described. Next, the Methods section describes how the pySUMMA can be used to create four different canopy shortwave radiation parameterizations of the Reynolds Mountain East catchment model. The Results section shows how to use pySUMMA and the Pandas library to reproduce Figure 1(above) from Clark et al. (2015).

Collectively, this Jupyter Notebook serves as an example of how hydrologic modeling can be conducted directly within a Jupyter Notebook by leveraging the pySUMMA library.

An important source of predictive differences among hydrologic and land-surface models is the method used to simulate the transmission and attenuation of shortwave radiation through the vegetation canopy. The main inter-model differences stem from (i) the methods used to simulate radiation transmission through homogenous vegetation [Dickinson, 1983; Sellers, 1985; Nijssen and Lettenmaier, 1999; Mahat and Tarboton, 2012]; (ii) the methods used to parameterize the impact of the canopy gap fraction on grid-average shortwave radiation fluxes [Cescatti, 1997; Kucharik et al., 1999; Niu and Yang, 2004; Essery et al., 2008]; and (iii) the methods used to represent spatial variability in vegetation type [Koster and Suarez, 1992; Bonan et al., 2002]. In this paper the parameterizations of canopy shortwave radiation are restricted to radiation transmission through homogenous vegetation, as this approach is used in many existing models. Recent advances in modeling the impact of canopy heterogeneity on grid average fluxes [e.g., Essery et al., 2008] are not included at this stage in model development, and will be considered in future work.

The methods considered for radiation transmission through homogenous vegetation allow for different levels of model complexity. At the simplest level we include an application of Beerâ€™s Law for direct-beam radiation (e.g., as used in VIC). At a more complex level, we include methods that model direct and diffuse beams separately and account for multiple scattering and multiple reflections [Nijssen and Lettenmaier, 1999]. Building additional complexity, we also include options for two-stream radiative transfer models as implemented in UEB [Mahat and Tarboton, 2012] and the Noah-MP model and CLM [Dickinson, 1983; Sellers, 1985; Oleson et al., 2010; Niu et al., 2011].

The interception of direct-beam shortwave radiation by vegetation at zenith angle ${\theta}_{zen}$ can be described as [e.g., Mahat and Tarboton, 2012]

\begin{equation*}
\frac{dQ_{swb}}{dz} = -Q_{swb}\frac{G\eta_{l}}{cos({\theta}_{zen})} \mspace{200mu} (1)
\end{equation*}

where ${dQ_{swb}}$ $(W {m}^{-2})$ is the downward shortwave radiation flux for the direct-beam, G (${m}^{2}{m}^{-2}$) is the leaf orientation factor defining the average area of leaves when viewed from direction ${\theta}_{zen}$, ${\eta_{l}}({m}^{2}{m}^{-3})$ is the leaf+stem density and z(m) is the vertical coordinate, positive downwards. Given the canopy depth ${D_{can}} = h_{top}^{veg} - h_{bot}^{veg}$ (m), integrating equation (1) results in Beerâ€™s law [Mahat and Tarboton, 2012]

\begin{equation*}
Q_{swb}(D_{can}) = Q_{swb}^{0} exp(-G\frac{\eta_{l} D_{can}}{cos({\theta}_{zen})}) \mspace{200mu} (2)
\end{equation*}

where ${dQ_{swb}^{0}}$ $(W {m}^{-2})$ defines the direct-beam shortwave radiation flux at the top of the canopy.

The penetration of direct beam shortwave radiation through the canopy in the absence of scattering can then be described as

\begin{equation*}
\tau_{pb} = \frac{Q_{swb}(D_{can})} {Q_{swb}^{0}} = exp(-G\frac{V_{ex}}{cos({\theta}_{zen})})\mspace{200mu} (3)
\end{equation*}

where $V_{ex} = {\eta_{l} D_{can}} ({m}^{2}{m}^{-2}) $ defines the total exposed leaf+stem area, that is, the leaf+stem area that is not buried by snow.

Under the simplifying assumption that the total transmission of direct and diffuse radiation (Ï„) is equal to the penetration of direct beam radiation, i.e.,

\begin{equation*}
\tau = \tau_{pb}\mspace{200mu} (4)
\end{equation*}

the shortwave radiation flux absorbed by the vegetation canopy and the ground surface is

\begin{equation*}
Q_{swnet}^{veg} = Q_{sw}^{0}(1-\tau)(1-\alpha^{veg})+Q_{sw}^{0} \tau \alpha^{sfc}(1-\tau)\mspace{200mu} (5)
\end{equation*}

\begin{equation*}
Q_{swnet}^{sfc} = Q_{sw}^{0}\tau(1-\alpha^{sfc})\mspace{200mu} (6)
\end{equation*}

where ${\alpha^{veg}}$ and $\alpha^{sfc}$(-) define the albedo of the vegetation canopy and ground surface respectively, and $Q_{sw}^{0} = Q_{swb}^{0} + Q_{swd}^{0} (W {m}^{-2})$ is the sum of direct-beam and diffuse radiation at the top of the canopy. Note that ${\alpha^{veg}} $(-) defines the bulk canopy albedo, which is typically much lower than the reflectance of an individual leaf due to the partial trapping of light by layers of leaves [Dickinson, 1983; Nijssen and Lettenmaier, 1999].

The first term on the right-hand-side of equation (5) represents the shortwave radiation absorbed by the canopy on the downward pass, while the second term represents the shortwave radiation that is reflected by the surface and absorbed on the upward pass. The approach described by equations (4) through (6) is similar to that used in the Variable Infiltration Capacity (VIC) model [Liang et al., 1994] and the Distributed Hydrology Soil Vegetation Model (DHSVM) [Wigmosta et al., 1994; Wigmosta and Lettenmaier, 1999], except that (i) VIC and DHSVM do not include a dependence on the solar zenith angle; and (ii) VIC assumes that all radiation reflected by the ground surface is lost through the top of the canopy, i.e., VIC does not include the second term in equation (5), which can be quite important when the ground surface is covered by snow.

The approach described above assumes all radiation is from a single direct beam and ignores the effects of scattering and multiple reflections. Nijssen and Lettenmaier [1999] showed that the penetration of diffuse radiation, $\tau_{pd}$, can be obtained by integrating over the upper hemisphere, as

\begin{equation*}
\tau_{pd} = \frac{1}{Q_{swd}^{0}} \int_{\Omega} {Q_{swd}^{0}}({\theta}_{zen}) \tau_{pb} ({\theta}_{zen}) cos({\theta}_{zen}) d\Omega \mspace{200mu} (7)
\end{equation*}

where ${Q_{swd}^{0}}({\theta}_{zen})$ is the above-canopy diffuse-beam radiation from the direction ${\theta}_{zen}$, and $\tau_{pb} ({\theta}_{zen})$ is the direct beam transmissivity from the direction ${\theta}_{zen}$ as defined in equation (16).

Assuming diffuse radiation is isotropic, a solution to equation (20) is [Nijssen and Lettenmaier, 1999; Mahat and Tarboton, 2012]

\begin{equation*}
\tau_{pd} = \frac{Q_{swb}(D_{can})} {Q_{swb}^{0}} = (1-GV_{ex})exp(-GV_{ex})+(GV_{ex})^2E_{i}(1,GV_{ex}) \mspace{200mu} (8)
\end{equation*}

where $E_{i}(n,x) $ is the exponential integral, and recall that $V_{ex}={\eta_{l} D_{can}}$ defines the total exposed leaf+stem area (${m}^{2}{m}^{-2}$), that is, the leaf+stem area that is not buried by snow. The total penetration of radiation through the canopy is then

\begin{equation*}
\tau_{pbd} = \frac{Q_{swb}^{0}}{Q_{sw}^{0}}\tau_{pb} + \frac{Q_{swd}^{0}}{Q_{sw}^{0}}\tau_{pd} \mspace{200mu} (9)
\end{equation*}

where $\tau_{pb} $ is given by the solution to Beerâ€™s law in equation (3).

Nijssen and Lettenmaier [1999] introduce a simple method to account for scattering and multiple reflections, providing total transmission $\tau$ as

\begin{equation*}
\tau = (\tau_{pbd})^{\varnothing_{s}} \mspace{200mu} (10)
\end{equation*}

where the exponent ${\varnothing_{s}}$ accounts for scattering of radiation within the vegetation canopy. The value of ${\varnothing_{s}}$ should generally fall in the range 0.7 to 0.85, and ${\varnothing_{s}}$ = 0.81 in the simulations presented by Nijssen and Lettenmaier [1999].

Given equation (10) and the assumption that all shortwave radiation reflected by the ground surface is diffuse, the shortwave radiation flux absorbed by the vegetation canopy and the ground surface can be given as

\begin{equation*}
Q_{swnet}^{veg} = Q_{sw}^{0}(1-\tau)(1-\alpha^{veg})+Q_{sw}^{0} m_{r} \tau \alpha^{sfc}(1-(\tau_{pd})^{\varnothing_{s}}) \mspace{200mu} (11)
\end{equation*}

\begin{equation*}
Q_{swnet}^{sfc} = Q_{sw}^{0} m_{r} \tau(1-\alpha^{sfc}) \mspace{200mu} (12)
\end{equation*}

where the terms on the right-hand-side of equation (11) describe the downward and upward radiation absorbed by the vegetation canopy, respectively, and $m_{r}$ accounts for multiple reflections between the ground and canopy.

Note that equation (11) assumes that all reflected radiation is diffuse, and reflection of upward radiation from the canopy is incorporated in the calculation of $\tau$ in equation (10) through the multiple reflection factor $m_{r}$.

The multiple reflection factor $m_{r}$ is derived from the infinite series of reflections between the surface and vegetation [Nijssen and Lettenmaier, 1999] as

\begin{equation*}
m_{r} = \sum_{i=0}^{\infty}[\alpha^{sfc}\alpha^{veg}(1-(\tau_{pd})^{\varnothing_{s}})]^{i} = \frac{1}{1-\alpha^{sfc}\alpha^{veg}(1-(\tau_{pd})^{\varnothing_{s}})} \mspace{200mu} (13)
\end{equation*}

with the factor $(1-(\tau_{pd})^{\varnothing_{s}})$ included to represent the fraction of radiation transmitted and scattered upwards through the canopy [see also Mahat and Tarboton, 2012].

Mahat and Tarboton [2012] describe a two stream radiative transfer model based on the assumptions that radiation is scattered equally in an upward and downward direction and that scattering is along the same path as the incoming radiation. These assumptions yield

\begin{equation*}
- \frac{dU_{sw}}{dz} = \frac{G\eta_{l}}{cos({\theta}_{zen})}[-U_{sw} + \alpha_{leaf} \frac{U_{sw}+Q_{sw}}{2}] \mspace{200mu} (14)
\end{equation*}

\begin{equation*}
- \frac{dQ_{sw}}{dz} = \frac{G\eta_{l}}{cos({\theta}_{zen})}[-Q_{sw} + \alpha_{leaf} \frac{U_{sw}+Q_{sw}}{2}] \mspace{200mu} (15)
\end{equation*}

where $\alpha_{leaf}$ (-) is the leaf scattering coefficient, distinguished from the bulk canopy albedo used in equation (5) and (11), and $U_{sw}$ and $Q_{sw}$ $(W {m}^{-2})$ are the intensity of the upward and downward beams. The solution for transmission over an infinitely deep canopy is [Mahat and Tarboton, 2012]

\begin{equation*}
\tau_{deep,b} = \frac{Q_{swb}(D_{can})} {Q_{swb}^{0}} = exp(-k^{'}G\frac{V_{ex}}{cos{\theta}_{zen}}) \mspace{200mu} (16)
\end{equation*}

\begin{equation*}
\tau_{deep,d} = \frac{Q_{swd}(D_{can})} {Q_{swd}^{0}} = (1-k^{'}GV_{ex})exp(-k^{'}GV_{ex})+(k^{'}GV_{ex})^2E_{i}(1,k^{'}GV_{ex}) \mspace{200mu} (17)
\end{equation*}

where $ k^{'} = \sqrt{1-\alpha_{leaf}} $ Equations (16) and (17) are similar to the expressions in equations (3) and (8), except the factor $ k^{'}$ is included to account for the effects of multiple scattering. Recall that $V_{ex}={\eta_{l}}D_{can}$.

The upward reflection factor $\alpha_{deep}^{veg} $ giving the fraction of radiation reflected from a deep canopy with multiple scattering is

\begin{equation*}
\alpha_{deep}^{veg} = \frac{U_{sw}(z)}{Q_{sw}(z)} = \frac{1-k^{'}}{1+k^{'}} \mspace{200mu} (18)
\end{equation*}

Mahat and Tarboton [2012] use the principle of superposition to derive solutions for a finite canopy, as

\begin{equation*}
\tau_{b} = \frac{Q_{swb}(D_{can})} {Q_{swb}^{0}} = \frac{\tau_{deep,b}[1-(\alpha_{deep}^{veg})^2]}{1-(\alpha_{deep}^{veg})^{2}(\tau_{deep,b})^2} \mspace{200mu} (19)
\end{equation*}

\begin{equation*}
\alpha_{b}^{veg} = \frac{U_{swb}^{0}} {Q_{swb}^{0}} = \frac{\alpha_{deep}^{veg}[1-(\tau_{deep,b})^2]}{1-(\alpha_{deep}^{veg})^{2}(\tau_{deep,b})^2} \mspace{200mu} (20)
\end{equation*}

where $\tau_{b}$ and $\alpha_{b}^{veg}$ define the transmittance and reflectance for direct-beam radiation. Equations (19) and (20) can also be used to obtain $\tau_{d}$ and $\alpha_{d}^{veg}$ by replacing the direct transmittance $\tau_{deep,b}$ with the diffuse transmittance $\tau_{deep,d}$ [Mahat and Tarboton, 2012].

The shortwave radiation flux absorbed by the vegetation canopy and the ground surface can now be given as (now distinguishing between direct and diffuse surface albedo, $\alpha_{b}^{sfc}$ and $\alpha_{d}^{sfc}$, as defined in the next section),

\begin{equation*}
Q_{swnet}^{veg} = [Q_{swb}^{0}(1-\tau_{b})(1-\alpha^{veg}_{b})+Q_{swd}^{0}(1-\tau_{d})(1-\alpha^{veg}_{d}] + (Q_{swb}^{0}\tau_{b}\alpha_{b}^{sfc}+Q_{swd}^{0}\tau_{d}\alpha_{d}^{sfc})m_{r}(1-\tau_{d}) \mspace{50mu} (21)
\end{equation*}

\begin{equation*}
Q_{swnet}^{sfc} = [Q_{swb}^{0}\tau_{b}(1-\alpha^{sfc}_{b})+Q_{swd}^{0}(\tau_{d})(1-\alpha^{sfc}_{d})]m_{r} \mspace{200mu} (22)
\end{equation*}

where $m_{r}$ is the multiple reflections factor as defined in equation (13) but computed using diffuse reflectances and defined as $m_{r} = [1-\alpha_{d}^{veg}\alpha_{d}^{sfc}(1-\tau_{d})]^{-1}$ . As in equation (11) it is assumed that all upward radiation is diffuse.

We also consider the two-stream radiative transfer model of Dickinson [1983] and Sellers [1985], as implemented in both the Community Land Model [Oleson et al., 2010] and the Noah-MP model [Niu et al., 2011]. In this approach fluxes are computed separately for visible and near-infra-red wavelengths. Complete algorithmic details are provided by Oleson et al. [2010] and are not repeated here.

The above description are taken from the Stomal Resistance Method section within the manual Structure for Unifying Multiple Modeling Alternatives (SUMMA), Version 1.0: Technical Description (April, 2015).

In [1]:

```
# User have to use "pysumma" Jupyter kernel, if you install pysumma conda virtual environment
import sys
!{sys.executable} -m pip install git+https://github.com/UW-Hydro/pysumma.git@master --upgrade
```

In [2]:

```
from pysumma import hydroshare_utils
import os
```

In [3]:

```
# Download SUMMA Model Instance from HydroShare
resource_id = '67d3a9d949cb4269b951f1c17cb751dc'
instance = hydroshare_utils.get_hs_resource(resource_id, os.getcwd())
```

In [4]:

```
!cd {instance}/; chmod +x ./installTestCases_local.sh
!cd {instance}/; ./installTestCases_local.sh
```

In [5]:

```
%pylab inline
import cartopy
import geoviews as gv
import holoviews as hv
hv.notebook_extension('bokeh')
```