add linke turbidity formulas module to atmosphere (#278)

* add linke turbidity formulas module

* add kasten 1996 pyrheliometric formula for linke turbidity factor
* move atmosphere into package, just for fun, why not?
* make __init__ with same api as before

Signed-off-by: Mark Mikofski <[email protected]>

* ignore venv, sublime

* revert api changes from 6e9aaf0

* append linke_turb_forms.py to atmosphere.py
* rm -rf atmosphere/

Signed-off-by: Mark Mikofski <[email protected]>

* move demo to tests

* import xrange from six for python 3 compatibility
* fix pyrheliometric formula is from 1980

* fix docstring

* just use `range` instead of importing `six`
* change function to `kasten80_lt` to match `gueymard94_pw`
* start adding references

* test kaste96_lt against linke turbidity using tmy3 at 700nm

* add lots of references
* cite text quoted directly from Ineichen
* acknowledge Armel
* fix docstrings for numpydocs

* removing filters from kasten pyrhelometric formula

* refactor kasten96_lt api simpler

* remove alpha argument
* update docstring to explain that if molineaux method, then aod is
shape is only aod at 700[nm] but if bird-hulstrom, then aod should be
both aod380 and aod500
* add see also to angstrom_aod_at_lambda
* add warning for air mass or pwat out of range
* check for bird-hulstrom, raise type error if neither molineaux nor
* also simplify angstrom_aod_at_lambda() api, change aod to aod0,
change lambda0 to lambda1, and add lambda0, update docstring
* update test_atmosphere.test_kasten96_lt with new changes

* add angstrom_alpha() method in atmosphere to calculate alpha

* simplify kasten96_lt API

* separate arguments so their shapes don't change depending on method,
ie: before aod could be 1, 2, N, N x 1 or N x 2 depending on whether
Molineaux or Bird-Hulstrom were used, now, there are 3 separate args,
 aod700, aod380 and aod500.
* fix misspelling of Hulstrom, no "d", note several sources also make
 this error :(
* add new methods to api.rst
* add test_angstrom tests
* add raises value error test for bad method
* add lower() to compare only lower case of method name

* add new linke and angstrom functions to api.rst

* fix more hulstrom errors, no "d"
* fix typos in referecnes

* add doi to sphinx docs as extlink

* add 'doi': ('http://dx/doi.org/%s', 'DOI: ') to extlink in conf.py
* change doi links everywhere to use :doi:`doi number` instead of links

* make landscape happy and use r"""docstring with :math:`\alpha`"""

* instead of """docstring with :math:`\\alpha`
* add Co-Action publishing for original Tellus A DOI before moving to
 Taylor & Francis

* reformat long lines to make landscape happy :)

* use broadband AOD arg instead of methods

* change kasten96_lt AOD args to just a single arg called aod_bb for
broadband AOD
* remove "method" arg and if-then-else on method.lower() in
kasten96_lt()
* add note that according to Molineaux, aod at 700 nm is a good
approximation for broadband aod
* add new function bird_hulstrom80_aod_bb to atmosphere.py
* add it to api.rst
* add new test for bird_hulstrom80_aod_bb() to test_atmosphere.py
* update test_kasten96_lt() to use aod_bb
* use verbose names:
  - s/am/airmass_absolute/g
  - s/pwat/precipitable_water/g
* fix angstrom aod test to use more realistic values
* don't use brackets around units
* remove test for raises ValueError in test_kasten

* use hard coded test values instead of TMY3

* values span pwat range 0 - 5 cm, and airmass between 1 & 2.5 atm
* remove tmy3 file
* don't try Bird-Hulstrom, not necessary
* make pylint happy, isclose(a, b) a & b should be array or iterable

* remove unused imports

* fix test kasten

* fix assert np.isclose() should be scalar only
* super simplify test_kasten96_lt() to ranges of inputs and expected
outputs, don't need to compare to actual data, just test that method
is working nominally

* remove clearsky import from test_atmosphere

* add missing assert statements for tests

* update hulstrom expected value

* remove old melbourne-fl test data

Author: mikofski

.gitignore

@@ -22,6 +22,7 @@ var/
 Drafts/
 */build/
 build/
+venv/
 
 *.egg-info/
 .installed.cfg
@@ -58,6 +59,8 @@ coverage.xml
 .pydevproject
 .spyderproject
 .idea/
+*.sublime-project
+*.sublime-workspace
 
 # Rope
 .ropeproject

docs/sphinx/source/api.rst

@@ -87,6 +87,10 @@ Airmass and atmospheric models
    atmosphere.alt2pres
    atmosphere.gueymard94_pw
    atmosphere.first_solar_spectral_correction
+   atmosphere.bird_hulstrom80_aod_bb
+   atmosphere.kasten96_lt
+   atmosphere.angstrom_aod_at_lambda
+   atmosphere.angstrom_alpha
 
 
 Irradiance

docs/sphinx/source/conf.py

@@ -270,7 +270,8 @@ def setup(app):
 extlinks = {'issue': ('https://github.com/pvlib/pvlib-python/issues/%s',
                       'GH'),
             'wiki': ('https://github.com/pvlib/pvlib-python/wiki/%s',
-                     'wiki ')}
+                     'wiki '),
+            'doi': ('http://dx.doi.org/%s', 'DOI: ')}
 
 # -- Options for manual page output ---------------------------------------
 

docs/sphinx/source/whatsnew/v0.4.4.txt

@@ -6,6 +6,9 @@ v0.4.4 (January xx, 2017)
 Enhancements
 ~~~~~~~~~~~~
 
+* Added Kasten pyrheliometric formula to calculate Linke turbidity factors with
+  improvements by Ineichen and Perez to extend range of air mass (:issue:`278`)
+
 
 Documentation
 ~~~~~~~~~~~~~
@@ -19,3 +22,4 @@ Contributors
 ~~~~~~~~~~~~
 
 * Will Holmgren
+* Mark Mikofski

pvlib/atmosphere.py

@@ -292,7 +292,7 @@ def gueymard94_pw(temp_air, relative_humidity):
     ----------
     .. [1] W. M. Keogh and A. W. Blakers, Accurate Measurement, Using Natural
        Sunlight, of Silicon Solar Cells, Prog. in Photovoltaics: Res.
-       and Appl. 2004, vol 12, pp. 1-19 (DOI: 10.1002/pip.517)
+       and Appl. 2004, vol 12, pp. 1-19 (:doi:`10.1002/pip.517`)
 
     .. [2] C. Gueymard, Analysis of Monthly Average Atmospheric Precipitable
        Water and Turbidity in Canada and Northern United States,
@@ -467,3 +467,193 @@ def first_solar_spectral_correction(pw, airmass_absolute, module_type=None,
         coeff[4]*np.sqrt(pw) + coeff[5]*ama/np.sqrt(pw))
 
     return modifier
+
+
+def bird_hulstrom80_aod_bb(aod380, aod500):
+    """
+    Approximate broadband aerosol optical depth.
+
+    Bird and Hulstrom developed a correlation for broadband aerosol optical
+    depth (AOD) using two wavelengths, 380 nm and 500 nm.
+
+    Parameters
+    ----------
+    aod380 : numeric
+        AOD measured at 380 nm
+    aod500 : numeric
+        AOD measured at 500 nm
+
+    Returns
+    -------
+    aod_bb : numeric
+        broadband AOD
+
+    See also
+    --------
+    kasten96_lt
+
+    References
+    ----------
+    [1] Bird and Hulstrom, "Direct Insolation Models" (1980)
+    `SERI/TR-335-344 <http://www.nrel.gov/docs/legosti/old/344.pdf>`_
+
+    [2] R. E. Bird and R. L. Hulstrom, "Review, Evaluation, and Improvement of
+    Direct Irradiance Models", Journal of Solar Energy Engineering 103(3),
+    pp. 182-192 (1981)
+    :doi:`10.1115/1.3266239`
+    """
+    # approximate broadband AOD using (Bird-Hulstrom 1980)
+    return 0.27583 * aod380 + 0.35 * aod500
+
+
+def kasten96_lt(airmass_absolute, precipitable_water, aod_bb):
+    """
+    Calculate Linke turbidity factor using Kasten pyrheliometric formula.
+
+    Note that broadband aerosol optical depth (AOD) can be approximated by AOD
+    measured at 700 nm according to Molineaux [4] . Bird and Hulstrom offer an
+    alternate approximation using AOD measured at 380 nm and 500 nm.
+
+    Based on original implementation by Armel Oumbe.
+
+    .. warning::
+        These calculations are only valid for air mass less than 5 atm and
+        precipitable water less than 5 cm.
+
+    Parameters
+    ----------
+    airmass_absolute : numeric
+        airmass, pressure corrected in atmospheres
+    precipitable_water : numeric
+        precipitable water or total column water vapor in centimeters
+    aod_bb : numeric
+        broadband AOD
+
+    Returns
+    -------
+    lt : numeric
+        Linke turbidity
+
+    See also
+    --------
+    bird_hulstrom80_aod_bb
+    angstrom_aod_at_lambda
+
+    References
+    ----------
+    [1] F. Linke, "Transmissions-Koeffizient und Trubungsfaktor", Beitrage
+    zur Physik der Atmosphare, Vol 10, pp. 91-103 (1922)
+
+    [2] F. Kasten, "A simple parameterization of the pyrheliometric formula for
+    determining the Linke turbidity factor", Meteorologische Rundschau 33,
+    pp. 124-127 (1980)
+
+    [3] Kasten, "The Linke turbidity factor based on improved values of the
+    integral Rayleigh optical thickness", Solar Energy, Vol. 56, No. 3,
+    pp. 239-244 (1996)
+    :doi:`10.1016/0038-092X(95)00114-7`
+
+    [4] B. Molineaux, P. Ineichen, N. O'Neill, "Equivalence of pyrheliometric
+    and monochromatic aerosol optical depths at a single key wavelength",
+    Applied Optics Vol. 37, issue 10, 7008-7018 (1998)
+    :doi:`10.1364/AO.37.007008`
+
+    [5] P. Ineichen, "Conversion function between the Linke turbidity and the
+    atmospheric water vapor and aerosol content", Solar Energy 82,
+    pp. 1095-1097 (2008)
+    :doi:`10.1016/j.solener.2008.04.010`
+
+    [6] P. Ineichen and R. Perez, "A new airmass independent formulation for
+    the Linke Turbidity coefficient", Solar Energy, Vol. 73, no. 3, pp. 151-157
+    (2002)
+    :doi:`10.1016/S0038-092X(02)00045-2`
+    """
+    # "From numerically integrated spectral simulations done with Modtran
+    # (Berk, 1989), Molineaux (1998) obtained for the broadband optical depth
+    # of a clean and dry atmospshere (fictitious atmosphere that comprises only
+    # the effects of Rayleigh scattering and absorption by the atmosphere gases
+    # other than the water vapor) the following expression"
+    # - P. Ineichen (2008)
+    delta_cda = -0.101 + 0.235 * airmass_absolute ** (-0.16)
+    # "and the broadband water vapor optical depth where pwat is the integrated
+    # precipitable water vapor content of the atmosphere expressed in cm and am
+    # the optical air mass. The precision of these fits is better than 1% when
+    # compared with Modtran simulations in the range 1 < am < 5 and
+    # 0 < pwat < 5 cm at sea level" - P. Ineichen (2008)
+    delta_w = 0.112 * airmass_absolute ** (-0.55) * precipitable_water ** 0.34
+    # broadband AOD
+    delta_a = aod_bb
+    # "Then using the Kasten pyrheliometric formula (1980, 1996), the Linke
+    # turbidity at am = 2 can be written. The extension of the Linke turbidity
+    # coefficient to other values of air mass was published by Ineichen and
+    # Perez (2002)" - P. Ineichen (2008)
+    lt = -(9.4 + 0.9 * airmass_absolute) * np.log(
+        np.exp(-airmass_absolute * (delta_cda + delta_w + delta_a))
+    ) / airmass_absolute
+    # filter out of extrapolated values
+    return lt
+
+
+def angstrom_aod_at_lambda(aod0, lambda0, alpha, lambda1=700.0):
+    r"""
+    Get AOD at specified wavelength using Angstrom turbidity model.
+
+    Parameters
+    ----------
+    aod0 : numeric
+        aerosol optical depth (AOD) measured at known wavelength
+    lambda0 : numeric
+        wavelength in nanometers corresponding to ``aod0``
+    alpha : numeric
+        Angstrom :math:`\alpha` exponent corresponding to ``aod0``
+    lambda1 : numeric
+        desired wavelength in nanometers, defaults to 700 nm
+
+    Returns
+    -------
+    aod1 : numeric
+        AOD at desired wavelength, ``lambda1``
+
+    See also
+    --------
+    angstrom_alpha
+
+    References
+    ----------
+    [1] Anders Angstrom, "On the Atmospheric Transmission of Sun Radiation and
+    On Dust in the Air", Geografiska Annaler Vol. 11, pp. 156-166 (1929) JSTOR
+    :doi:`10.2307/519399`
+
+    [2] Anders Angstrom, "Techniques of Determining the Turbidity of the
+    Atmosphere", Tellus 13:2, pp. 214-223 (1961) Taylor & Francis
+    :doi:`10.3402/tellusa.v13i2.9493` and Co-Action Publishing
+    :doi:`10.1111/j.2153-3490.1961.tb00078.x`
+    """
+    return aod0 * ((lambda1 / lambda0) ** (-alpha))
+
+
+def angstrom_alpha(aod1, lambda1, aod2, lambda2):
+    r"""
+    Calculate Angstrom alpha exponent.
+
+    Parameters
+    ----------
+    aod1 : numeric
+        first aerosol optical depth
+    lambda1 : numeric
+        wavelength in nanometers corresponding to ``aod1``
+    aod2 : numeric
+        second aerosol optical depth
+    lambda2 : numeric
+        wavelength in nanometers corresponding to ``aod2``
+
+    Returns
+    -------
+    alpha : numeric
+        Angstrom :math:`\alpha` exponent for AOD in ``(lambda1, lambda2)``
+
+    See also
+    --------
+    angstrom_aod_at_lambda
+    """
+    return - np.log(aod1 / aod2) / np.log(lambda1 / lambda2)

pvlib/test/test_atmosphere.py

@@ -1,4 +1,3 @@
-import datetime
 import itertools
 
 import numpy as np
@@ -7,10 +6,8 @@
 import pytest
 from numpy.testing import assert_allclose
 
-from pvlib.location import Location
-from pvlib import solarposition
 from pvlib import atmosphere
-
+from pvlib import solarposition
 
 latitude, longitude, tz, altitude = 32.2, -111, 'US/Arizona', 700
 
@@ -113,3 +110,43 @@ def test_first_solar_spectral_correction_supplied():
 def test_first_solar_spectral_correction_ambiguous():
     with pytest.raises(TypeError):
         atmosphere.first_solar_spectral_correction(1, 1)
+
+
+def test_kasten96_lt():
+    """Test Linke turbidity factor calculated from AOD, Pwat and AM"""
+    amp = np.array([1, 3, 5])
+    pwat = np.array([0, 2.5, 5])
+    aod_bb = np.array([0, 0.1, 1])
+    lt_expected = np.array(
+        [[[1.3802, 2.4102, 11.6802],
+          [1.16303976, 2.37303976, 13.26303976],
+          [1.12101907, 2.51101907, 15.02101907]],
+
+         [[2.95546945, 3.98546945, 13.25546945],
+          [2.17435443, 3.38435443, 14.27435443],
+          [1.99821967, 3.38821967, 15.89821967]],
+
+         [[3.37410769, 4.40410769, 13.67410769],
+          [2.44311797, 3.65311797, 14.54311797],
+          [2.23134152, 3.62134152, 16.13134152]]]
+    )
+    lt = atmosphere.kasten96_lt(*np.meshgrid(amp, pwat, aod_bb))
+    assert np.allclose(lt, lt_expected, 1e-3)
+    return lt
+
+
+def test_angstrom_aod():
+    """Test Angstrom turbidity model functions."""
+    aod550 = 0.15
+    aod1240 = 0.05
+    alpha = atmosphere.angstrom_alpha(aod550, 550, aod1240, 1240)
+    assert np.isclose(alpha, 1.3513924317859232)
+    aod700 = atmosphere.angstrom_aod_at_lambda(aod550, 550, alpha)
+    assert np.isclose(aod700, 0.10828110997681031)
+
+
+def test_bird_hulstrom80_aod_bb():
+    """Test Bird_Hulstrom broadband AOD."""
+    aod380, aod500 = 0.22072480948195175, 0.1614279181106312
+    bird_hulstrom = atmosphere.bird_hulstrom80_aod_bb(aod380, aod500)
+    assert np.isclose(0.11738229553812768, bird_hulstrom)