horizonAdjustment.ipynb

horizonCorrection01.ipynb

{
  "cells": [
    {
      "metadata": {
        "toc": true
      },
      "cell_type": "markdown",
      "source": "<h1>Table of Contents<span class=\"tocSkip\"></span></h1>\n<div class=\"toc\"><ul class=\"toc-item\"><li><span><a href=\"#Baseline-$\\alpha$-approach\" data-toc-modified-id=\"Baseline-$\\alpha$-approach-1\"><span class=\"toc-item-num\">1&nbsp;&nbsp;</span>Baseline $\\alpha$ approach</a></span></li><li><span><a href=\"#Mean-of-the-set-of-$c$'s\" data-toc-modified-id=\"Mean-of-the-set-of-$c$'s-2\"><span class=\"toc-item-num\">2&nbsp;&nbsp;</span>Mean of the set of $c$'s</a></span></li><li><span><a href=\"#Direct-radiation:-horizon-shading\" data-toc-modified-id=\"Direct-radiation:-horizon-shading-3\"><span class=\"toc-item-num\">3&nbsp;&nbsp;</span>Direct radiation: horizon shading</a></span><ul class=\"toc-item\"><li><span><a href=\"#Visibility\" data-toc-modified-id=\"Visibility-3.1\"><span class=\"toc-item-num\">3.1&nbsp;&nbsp;</span>Visibility</a></span></li><li><span><a href=\"#Shading\" data-toc-modified-id=\"Shading-3.2\"><span class=\"toc-item-num\">3.2&nbsp;&nbsp;</span>Shading</a></span></li><li><span><a href=\"#Ratio\" data-toc-modified-id=\"Ratio-3.3\"><span class=\"toc-item-num\">3.3&nbsp;&nbsp;</span>Ratio</a></span></li><li><span><a href=\"#Double-series\" data-toc-modified-id=\"Double-series-3.4\"><span class=\"toc-item-num\">3.4&nbsp;&nbsp;</span>Double series</a></span></li></ul></li><li><span><a href=\"#Diffuse-radiation:-sky-view-factor\" data-toc-modified-id=\"Diffuse-radiation:-sky-view-factor-4\"><span class=\"toc-item-num\">4&nbsp;&nbsp;</span>Diffuse radiation: sky-view factor</a></span></li><li><span><a href=\"#Horizon-adjusted-$\\alpha$\" data-toc-modified-id=\"Horizon-adjusted-$\\alpha$-5\"><span class=\"toc-item-num\">5&nbsp;&nbsp;</span>Horizon-adjusted $\\alpha$</a></span></li></ul></div>"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "# Baseline $\\alpha$ approach\n[Bair et al 2015](https://people.eri.ucsb.edu/~nbair/storage/bair_et_al_2015_frontiers.pdf) give broadband snow albedo $\\alpha$ as:"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "$\\alpha = \\frac{D_\\uparrow}{cB_\\downarrow + D_\\downarrow}$"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "Where $B_\\downarrow$ and $D_\\downarrow$ are the direct and diffuse downwelling global radiation, and the terrain correction factor $c$ is calculated via $\\theta$ (the local illumination angle) and $\\theta_0$ (the illumination angle on a horizontal surface):"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "$c = \\frac{cos\\theta}{cos\\theta_0} $"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "# Mean of the set of $c$'s"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "In our faceted model, an overall planar fit is replaced by a large number of data subsets, each with their own planar fit.\n\nEach of these facets therefore has its own $c$-value. The mean of these values can be expressed as follows, with $n$ equal to the number of facets:"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "$\\overline c =\\frac{1}{n} \\displaystyle\\sum_{i=1}^{n} \\frac{cos\\theta_i}{cos\\theta_0} $"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "# Direct radiation: horizon shading\n"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "## Visibility"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "Regarding the 'visibility' of a given point $i$ from the horizon of another point $H_f(j)$, [Dozier et al 1981](http://www2.bren.ucsb.edu/~dozier/Pubs/DozierEtAlCompGeosci1981.pdf) point out:\n\n    \"A more efficient algorithm can be obtained by noting\n    that, in testing whether j is the horizon of point i, if\n    slope(i, j) > slope(j, Hf(j)) then all of the points forward\n    of j need not be checked because clearly they are not\n    'visible' from point i, and Hf(i)=j.\" (pg. 146-147)\n    \nHere $H_f(j)$ indicates the horizon of $j$ in the 'forward' $f$ direction, but this $f$ can be subsituted for any angle (i.e. azimuth) of our choosing."
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "## Shading"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "Based on the above, to determine whether a given point $j$ is shaded from direct solar radiation at a given time $t$, we can evaluate the truth of:\n\n$slope\\left(j, H_\\phi(j)\\right) > \\theta_t$\n\nHere, $H_\\phi$ is the horizon of $j$ in the direction of the solar azimuth $\\phi$, and $\\theta_t$ is the solar altitude at time $t$. \n\nFollowing Dozier et al, if the slope between $j$ and $H_\\phi(j)$ is greater than $\\theta_t$, then the point $j$ is shaded from direct solar radiation at time $t$."
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "## Ratio"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "For any real number $x$, the sign function: \n\n\n$\\Large \\frac{|x|}{x}$$ = \\begin{cases} -1 & \\text{for } x< 1\\\\ 0 & \\text{for } x=0\\\\ 1 & \\text{for } x> 1\\\\ \\end{cases}$"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "Therefore:\n\n$\\Large\n\\frac {|\\theta_t \\ - \\ slope\\left(j, H_\\phi(j)\\right)|}{\\theta_t \\ - \\ slope\\left(j, H_\\phi(j)\\right)} $$= \\begin{cases} -1 & \\text{if the point } j \\text{ is shaded from direct solar radiation by its horizon} \\\\ 0 & \\text{if slope between } j \\text{ and its horizon equals } \\theta_t \\\\ 1 & \\text{if the point } j \\text{ is exposed to direct solar radiation} \\\\ \\end{cases}\n$"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "For a range $n$ equal to the number of points on a given facet, the following sum returns the ratio* of the given facet which is visible to the sun at time $t$:"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "$\n\\frac{1}{2n} \\displaystyle\\sum_{j=1}^{n} \\frac {|\\theta_t \\ - \\ slope\\left(j, H_\\phi(i)\\right)|}{\\theta_t \\ - \\ slope\\left(j, H_\\phi(j) \\right)} + 1\n$"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "*Further thought can be given to how to account for the edge-case where $slope\\left(j, H_\\phi(i)\\right) = \\theta_t$."
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "## Double series"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "We recall that the mean value of $c$ across $n$ facets is:"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "$\\overline c =\\frac{1}{n} \\displaystyle\\sum_{i=1}^{n} \\frac{cos\\theta_i}{cos\\theta_0} $"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "We can further weight these terrain correction factors, to create a coefficient $c_H$, the 'horizon-adjusted terrain correction factor', by multiplying $c$ by the ratio of the given facet which is visible to the sun at time $t$. The mean of $c_H$ across a set of $m$ facets, each comprised of $n$ points, can be expressed as a [double series](https://mathworld.wolfram.com/DoubleSeries.html):"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "$\\overline {c_H} \n=\\left(\\frac{1}{m} \\displaystyle\\sum_{i=1}^{m} \\frac{cos\\theta_i}{cos\\theta_0}\\right) \\left( \\frac{1}{2n} \\displaystyle\\sum_{j=1}^{n} \\frac {|\\theta_t \\ - \\ slope\\left(j, H_\\phi(j)\\right)|}{\\theta_t \\ - \n\\ slope\\left(j, H_\\phi(j)\\right)} + 1 \\right) \n$"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "$\n\\ \\quad = \\frac{1}{2mn} \\left(\\displaystyle\\sum_{i=1}^{m} \\frac{cos\\theta_i}{cos\\theta_0}\\right) \\left( \\displaystyle\\sum_{j=1}^{n} \\frac {|\\theta_t - slope\\left(j, H_\\phi(j)\\right)|}{\\theta_t - slope\\left(j, H_\\phi(j)\\right)} + 1 \\right)\n$"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "This coefficient $\\overline{c_H}$ can be used in place of $c$ in the albedo equation, to incorporate horizon adjustment into the scaling of $B_\\downarrow$.\n\nIt is preferable to calculate $\\overline{c_H}$ as a double (i.e. nested) series, rather than calculating $\\overline c$ separately and then multiplying it by a second coefficient for shade ratio. This approach allows us to scale the specific values of $c$ which are shaded, rather than applying a blanket shade ratio to all $c$'s.\n\n(There are a few different options for dealing with facets for which $c_H$ cannot be calculated; e.g. facets obstructed by the radiometer boom. One choice would be to fill those values with some sort of default or mean value. Another would be to exclude these facets from the calculation of $\\overline{c_H}$ and simply take the number of facets, $m$, to be the number of facets for which $c_H$ is calcuable.)"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "# Diffuse radiation: sky-view factor\n\n[Dozier and Frew 1979](http://www2.bren.ucsb.edu/~dozier/Pubs/DozierFrewIEEE1990.pdf) Eq. 7b provides an approximation of the the sky-view factor $V_d$ for a given point on a terrain grid:"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "$V_d \\approx \\frac{1}{2\\pi} \\int_{0}^{2\\pi} \\left[\\cos S \\ \\sin^2 H_\\phi + \\sin S \\ \\cos (\\phi - A) \\ \\times (H_\\phi - \\sin H_\\phi \\cos H_\\phi)\\right] d\\phi$"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "The mean sky-view factor $\\overline{V_d}$ for a terrain grid of $n$ points can be expressed as follows:"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "$\\overline{V_d} \\approx \\frac{1}{n} \\displaystyle\\sum_{i=1}^{n} \\frac{1}{2\\pi} \\int_{0}^{2\\pi} \\left[\\cos S \\ \\sin^2 H_\\phi + \\sin S \\ \\cos (\\phi - A) \\ \\times (H_\\phi - \\sin H_\\phi \\cos H_\\phi)\\right] d\\phi$"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "This coefficient $\\overline{V_d}$ can be used to scale diffuse downwelling radiation $D_\\downarrow$ in the albedo equation."
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "# Horizon-adjusted $\\alpha$"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "We recall that the baseline formula for broadband snow albedo $\\alpha$ is:"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "$\\alpha = \\frac{D_\\uparrow}{c B_\\downarrow + D_\\downarrow}$"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "To account for horizon shading of direct radiation, we can substitute $\\overline{c_H}$, defined above, for $c$. We can further multiply $D_\\downarrow$ by a mean sky-view factor $\\overline{V_d}$ to appropriately scale incoming diffuse radiation. This results in a new equation for $\\alpha$, adjusted for horizon shading and sky-view factors:"
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "$\\Large \\alpha = \\frac{D_\\uparrow}{\\overline {c_H} B_\\downarrow + \\overline {V_d}D_\\downarrow} $"
    }
  ],
  "metadata": {
    "kernelspec": {
      "name": "lidar20200821",
      "display_name": "LIDAR20200821",
      "language": "python"
    },
    "language_info": {
      "name": "python",
      "version": "3.7.7",
      "mimetype": "text/x-python",
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "pygments_lexer": "ipython3",
      "nbconvert_exporter": "python",
      "file_extension": ".py"
    },
    "toc": {
      "nav_menu": {},
      "number_sections": true,
      "sideBar": true,
      "skip_h1_title": false,
      "base_numbering": 1,
      "title_cell": "Table of Contents",
      "title_sidebar": "Contents",
      "toc_cell": true,
      "toc_position": {},
      "toc_section_display": true,
      "toc_window_display": true
    },
    "varInspector": {
      "window_display": false,
      "cols": {
        "lenName": 16,
        "lenType": 16,
        "lenVar": 40
      },
      "kernels_config": {
        "python": {
          "library": "var_list.py",
          "delete_cmd_prefix": "del ",
          "delete_cmd_postfix": "",
          "varRefreshCmd": "print(var_dic_list())"
        },
        "r": {
          "library": "var_list.r",
          "delete_cmd_prefix": "rm(",
          "delete_cmd_postfix": ") ",
          "varRefreshCmd": "cat(var_dic_list()) "
        }
      },
      "types_to_exclude": [
        "module",
        "function",
        "builtin_function_or_method",
        "instance",
        "_Feature"
      ]
    },
    "gist": {
      "id": "c60bef384331a8e4bad19baff2a31bd9",
      "data": {
        "description": "horizonAdjustment.ipynb",
        "public": true
      }
    },
    "_draft": {
      "nbviewer_url": "https://gist.github.com/c60bef384331a8e4bad19baff2a31bd9"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 4
}