diff --git a/.travis.yml b/.travis.yml index 944cbb8e21..7da069cf00 100644 --- a/.travis.yml +++ b/.travis.yml @@ -23,8 +23,8 @@ addons: - ccache cache: - directories: - - $HOME/.ccache +# - pip + - ccache # setup miniconda for numpy, scipy, pandas before_install: @@ -46,6 +46,14 @@ install: - echo "install" - conda env create --file ci/requirements-$CONDA_ENV.yml - source activate test_env # all envs are named test_env in the yml files + # needed to make sure that pandas is compiled against the right + # version of numpy + - if [[ "$CONDA_ENV" == "py27-min" ]]; then + pip uninstall numpy --yes; + pip uninstall pandas --yes; + pip install --no-cache-dir numpy==1.9.0; + pip install --no-cache-dir pandas==0.14.0; + fi - conda list - echo $PATH - ls -l /home/travis/miniconda/envs/test_env/lib diff --git a/ci/requirements-py27-min.yml b/ci/requirements-py27-min.yml index 5d24dc4d83..d132091fdc 100644 --- a/ci/requirements-py27-min.yml +++ b/ci/requirements-py27-min.yml @@ -1,8 +1,6 @@ name: test_env dependencies: - python=2.7 - - numpy==1.8.2 - - pandas==0.13.1 - pytz - pytest - pytest-cov diff --git a/docs/sphinx/source/whatsnew/v0.4.0.txt b/docs/sphinx/source/whatsnew/v0.4.0.txt index d5adb5853c..7a43bccfc6 100644 --- a/docs/sphinx/source/whatsnew/v0.4.0.txt +++ b/docs/sphinx/source/whatsnew/v0.4.0.txt @@ -14,6 +14,10 @@ API Changes * Remove unneeded module argument from singlediode function. (:issue:`200`) * In ``pvlib.irradiance.perez``, renamed argument ``modelt`` to ``model``. (:issue:`196`) +* Functions in the irradiance module now work with scalar inputs + in addition to arrays and Series. Furthermore, these functions no + longer promote scalar or array input to Series output. + Also applies to atmosphere.relativeairmass. (:issue:`201`, :issue:`214`) Enhancements @@ -22,6 +26,8 @@ Enhancements * Adds the First Solar spectral correction model. (:issue:`115`) * Adds the Gueymard 1994 integrated precipitable water model. (:issue:`115`) * Adds the PVWatts DC, AC, and system losses model. (:issue:`195`) +* Improve PEP8 conformity in irradiance module. (:issue:`214`) +* irradiance.disc is up to 10x faster. (:issue:`214`) Bug fixes @@ -44,6 +50,17 @@ Other * Switch to the py.test testing framework. (:issue:`204`) * Reconfigure Appveyor CI builds and resolve an issue in which the command line length was too long. (:issue:`207`) +* Manually build numpy and pandas for the min requirements test. + This is needed to avoid Continuum's bad practice of bundling scipy + with pandas. (:issue:`214`) + + +Requirements +~~~~~~~~~~~~ + +* pvlib now requires pandas >= 0.14.0 and numpy >= 1.9.0, both + released in 2014. Most of pvlib will work with lesser versions. + (:issue:`214`) Code Contributors diff --git a/pvlib/atmosphere.py b/pvlib/atmosphere.py index 494e0dcc22..d5e8cbc1d3 100644 --- a/pvlib/atmosphere.py +++ b/pvlib/atmosphere.py @@ -6,6 +6,7 @@ from __future__ import division import numpy as np +import pandas as pd from warnings import warn APPARENT_ZENITH_MODELS = ('simple', 'kasten1966', 'kastenyoung1989', @@ -20,12 +21,12 @@ def pres2alt(pressure): Parameters ---------- - pressure : scalar or Series + pressure : numeric Atmospheric pressure (Pascals) Returns ------- - altitude : scalar or Series + altitude : numeric Altitude in meters above sea level Notes @@ -45,12 +46,12 @@ def pres2alt(pressure): References ----------- - - "A Quick Derivation relating altitude to air pressure" from Portland - State Aerospace Society, Version 1.03, 12/22/2004. + [1] "A Quick Derivation relating altitude to air pressure" from + Portland State Aerospace Society, Version 1.03, 12/22/2004. ''' alt = 44331.5 - 4946.62 * pressure ** (0.190263) + return alt @@ -60,12 +61,12 @@ def alt2pres(altitude): Parameters ---------- - Altitude : scalar or Series + altitude : numeric Altitude in meters above sea level Returns ------- - Pressure : scalar or Series + pressure : numeric Atmospheric pressure (Pascals) Notes @@ -85,9 +86,8 @@ def alt2pres(altitude): References ----------- - - "A Quick Derivation relating altitude to air pressure" from Portland - State Aerospace Society, Version 1.03, 12/22/2004. + [1] "A Quick Derivation relating altitude to air pressure" from + Portland State Aerospace Society, Version 1.03, 12/22/2004. ''' press = 100 * ((44331.514 - altitude) / 11880.516) ** (1 / 0.1902632) @@ -111,16 +111,15 @@ def absoluteairmass(airmass_relative, pressure=101325.): Parameters ---------- - - airmass_relative : scalar or Series + airmass_relative : numeric The airmass at sea-level. - pressure : scalar or Series + pressure : numeric The site pressure in Pascal. Returns ------- - scalar or Series + airmass_absolute : numeric Absolute (pressure corrected) airmass References @@ -128,7 +127,6 @@ def absoluteairmass(airmass_relative, pressure=101325.): [1] C. Gueymard, "Critical analysis and performance assessment of clear sky solar irradiance models using theoretical and measured data," Solar Energy, vol. 51, pp. 121-138, 1993. - ''' airmass_absolute = airmass_relative * pressure / 101325. @@ -147,15 +145,14 @@ def relativeairmass(zenith, model='kastenyoung1989'): Parameters ---------- - - zenith : float or Series + zenith : numeric Zenith angle of the sun in degrees. Note that some models use the apparent (refraction corrected) zenith angle, and some models use the true (not refraction-corrected) zenith angle. See model descriptions to determine which type of zenith angle is required. Apparent zenith angles must be calculated at sea level. - model : String + model : string Available models include the following: * 'simple' - secant(apparent zenith angle) - @@ -175,13 +172,12 @@ def relativeairmass(zenith, model='kastenyoung1989'): Returns ------- - airmass_relative : float or Series - Relative airmass at sea level. Will return NaN values for any + airmass_relative : numeric + Relative airmass at sea level. Will return NaN values for any zenith angle greater than 90 degrees. References ---------- - [1] Fritz Kasten. "A New Table and Approximation Formula for the Relative Optical Air Mass". Technical Report 136, Hanover, N.H.: U.S. Army Material Command, CRREL. @@ -207,7 +203,9 @@ def relativeairmass(zenith, model='kastenyoung1989'): Sandia Report, (2012). ''' - z = zenith + # need to filter first because python 2.7 does not support raising a + # negative number to a negative power. + z = np.where(zenith > 90, np.nan, zenith) zenith_rad = np.radians(z) model = model.lower() @@ -237,10 +235,8 @@ def relativeairmass(zenith, model='kastenyoung1989'): else: raise ValueError('%s is not a valid model for relativeairmass', model) - try: - am[z > 90] = np.nan - except TypeError: - am = np.nan if z > 90 else am + if isinstance(zenith, pd.Series): + am = pd.Series(am, index=zenith.index) return am @@ -282,14 +278,14 @@ def gueymard94_pw(temp_air, relative_humidity): Parameters ---------- - temp_air : array-like + temp_air : numeric ambient air temperature at the surface (C) - relative_humidity : array-like + relative_humidity : numeric relative humidity at the surface (%) Returns ------- - pw : array-like + pw : numeric precipitable water (cm) References @@ -349,7 +345,7 @@ def first_solar_spectral_correction(pw, airmass_absolute, module_type=None, Pwat where: * 0.5 cm <= Pwat <= 5 cm - * 1.0 <= AMa <= 5.0 + * 1.0 <= AMa <= 5.0 * Spectral range is limited to that of CMP11 (280 nm to 2800 nm) * spectrum simulated on a plane normal to the sun * All other parameters fixed at G173 standard @@ -408,8 +404,8 @@ def first_solar_spectral_correction(pw, airmass_absolute, module_type=None, radiative transfer of sunshine: algorithms and performance assessment. Cocoa, FL: Florida Solar Energy Center, 1995. .. [2] Lee, Mitchell, and Panchula, Alex. "Spectral Correction for - Photovoltaic Module Performance Based on Air Mass and Precipitable - Water." IEEE Photovoltaic Specialists Conference, Portland, 2016 + Photovoltaic Module Performance Based on Air Mass and Precipitable + Water." IEEE Photovoltaic Specialists Conference, Portland, 2016 .. [3] Marion, William F., et al. User's Manual for Data for Validating Models for PV Module Performance. National Renewable Energy Laboratory, 2014. http://www.nrel.gov/docs/fy14osti/61610.pdf @@ -423,30 +419,29 @@ def first_solar_spectral_correction(pw, airmass_absolute, module_type=None, if np.min(pw) < 0.1: pw = np.maximum(pw, 0.1) - warn('Exceptionally low Pwat values replaced with 0.1 cm to prevent'+ - ' model divergence') + warn('Exceptionally low Pwat values replaced with 0.1 cm to prevent' + + ' model divergence') # Warn user about Pwat data that is exceptionally high if np.max(pw) > 8: warn('Exceptionally high Pwat values. Check input data:' + - ' model may diverge in this range') + ' model may diverge in this range') # *** AMa *** # Replace Extremely High AM with AM 10 to prevent model divergence # AM > 10 will only occur very close to sunset - if np.max(airmass_absolute) > 10: - airmass_absolute = np.minimum(airmass_absolute,10) - + if np.max(airmass_absolute) > 10: + airmass_absolute = np.minimum(airmass_absolute, 10) + # Warn user about AMa data that is exceptionally low if np.min(airmass_absolute) < 0.58: warn('Exceptionally low air mass: ' + - 'model not intended for extra-terrestrial use') - # pvl_absoluteairmass(1,pvl_alt2pres(4340)) = 0.58 - # Elevation of Mina Pirquita, Argentian = 4340 m. Highest elevation city - # with population over 50,000. - + 'model not intended for extra-terrestrial use') + # pvl_absoluteairmass(1,pvl_alt2pres(4340)) = 0.58 Elevation of + # Mina Pirquita, Argentian = 4340 m. Highest elevation city with + # population over 50,000. _coefficients = {} _coefficients['cdte'] = ( @@ -468,9 +463,9 @@ def first_solar_spectral_correction(pw, airmass_absolute, module_type=None, # Evaluate Spectral Shift coeff = coefficients - AMa = airmass_absolute + ama = airmass_absolute modifier = ( - coeff[0] + coeff[1]*AMa + coeff[2]*pw + coeff[3]*np.sqrt(AMa) + - + coeff[4]*np.sqrt(pw) + coeff[5]*AMa/np.sqrt(pw)) + coeff[0] + coeff[1]*ama + coeff[2]*pw + coeff[3]*np.sqrt(ama) + + coeff[4]*np.sqrt(pw) + coeff[5]*ama/np.sqrt(pw)) return modifier diff --git a/pvlib/irradiance.py b/pvlib/irradiance.py index fd4cf3ee8e..bb27f4faa7 100644 --- a/pvlib/irradiance.py +++ b/pvlib/irradiance.py @@ -1,22 +1,24 @@ """ -The ``irradiance`` module contains functions for modeling -global horizontal irradiance, direct normal irradiance, -diffuse horizontal irradiance, and total irradiance -under various conditions. +The ``irradiance`` module contains functions for modeling global +horizontal irradiance, direct normal irradiance, diffuse horizontal +irradiance, and total irradiance under various conditions. """ from __future__ import division import logging -pvl_logger = logging.getLogger('pvlib') import datetime +from collections import OrderedDict import numpy as np import pandas as pd from pvlib import tools from pvlib import solarposition +from pvlib import atmosphere + +pvl_logger = logging.getLogger('pvlib') SURFACE_ALBEDOS = {'urban': 0.18, 'grass': 0.20, @@ -32,9 +34,6 @@ 'fresh steel': 0.35, 'dirty steel': 0.08} -# would be nice if this took pandas index as well. -# Use try:except of isinstance. - def extraradiation(datetime_or_doy, solar_constant=1366.1, method='spencer'): """ @@ -174,30 +173,32 @@ def _doy_to_timestamp(doy, epoch='2013-12-31'): def aoi_projection(surface_tilt, surface_azimuth, solar_zenith, solar_azimuth): """ - Calculates the dot product of the solar vector and the surface normal. + Calculates the dot product of the solar vector and the surface + normal. Input all angles in degrees. Parameters - ========== - - surface_tilt : float or Series. + ---------- + surface_tilt : numeric Panel tilt from horizontal. - surface_azimuth : float or Series. + surface_azimuth : numeric Panel azimuth from north. - solar_zenith : float or Series. + solar_zenith : numeric Solar zenith angle. - solar_azimuth : float or Series. + solar_azimuth : numeric Solar azimuth angle. Returns - ======= - float or Series. Dot product of panel normal and solar angle. + ------- + projection : numeric + Dot product of panel normal and solar angle. """ projection = ( - tools.cosd(surface_tilt) * tools.cosd(solar_zenith) + tools.sind(surface_tilt) * - tools.sind(solar_zenith) * tools.cosd(solar_azimuth - surface_azimuth)) + tools.cosd(surface_tilt) * tools.cosd(solar_zenith) + + tools.sind(surface_tilt) * tools.sind(solar_zenith) * + tools.cosd(solar_azimuth - surface_azimuth)) try: projection.name = 'aoi_projection' @@ -215,23 +216,24 @@ def aoi(surface_tilt, surface_azimuth, solar_zenith, solar_azimuth): Input all angles in degrees. Parameters - ========== - - surface_tilt : float or Series. + ---------- + surface_tilt : numeric Panel tilt from horizontal. - surface_azimuth : float or Series. + surface_azimuth : numeric Panel azimuth from north. - solar_zenith : float or Series. + solar_zenith : numeric Solar zenith angle. - solar_azimuth : float or Series. + solar_azimuth : numeric Solar azimuth angle. Returns - ======= - float or Series. Angle of incidence in degrees. + ------- + aoi : numeric + Angle of incidence in degrees. """ - projection = aoi_projection(surface_tilt, surface_azimuth, solar_zenith, solar_azimuth) + projection = aoi_projection(surface_tilt, surface_azimuth, + solar_zenith, solar_azimuth) aoi_value = np.rad2deg(np.arccos(projection)) try: @@ -245,27 +247,27 @@ def aoi(surface_tilt, surface_azimuth, solar_zenith, solar_azimuth): def poa_horizontal_ratio(surface_tilt, surface_azimuth, solar_zenith, solar_azimuth): """ - Calculates the ratio of the beam components of the - plane of array irradiance and the horizontal irradiance. + Calculates the ratio of the beam components of the plane of array + irradiance and the horizontal irradiance. Input all angles in degrees. Parameters - ========== - - surface_tilt : float or Series. + ---------- + surface_tilt : numeric Panel tilt from horizontal. - surface_azimuth : float or Series. + surface_azimuth : numeric Panel azimuth from north. - solar_zenith : float or Series. + solar_zenith : numeric Solar zenith angle. - solar_azimuth : float or Series. + solar_azimuth : numeric Solar azimuth angle. Returns - ======= - float or Series. Ratio of the plane of array irradiance to the - horizontal plane irradiance + ------- + ratio : numeric + Ratio of the plane of array irradiance to the horizontal plane + irradiance """ cos_poa_zen = aoi_projection(surface_tilt, surface_azimuth, @@ -284,31 +286,32 @@ def poa_horizontal_ratio(surface_tilt, surface_azimuth, return ratio -def beam_component(surface_tilt, surface_azimuth, - solar_zenith, solar_azimuth, dni): +def beam_component(surface_tilt, surface_azimuth, solar_zenith, solar_azimuth, + dni): """ Calculates the beam component of the plane of array irradiance. Parameters ---------- - surface_tilt : float or Series. + surface_tilt : numeric Panel tilt from horizontal. - surface_azimuth : float or Series. + surface_azimuth : numeric Panel azimuth from north. - solar_zenith : float or Series. + solar_zenith : numeric Solar zenith angle. - solar_azimuth : float or Series. + solar_azimuth : numeric Solar azimuth angle. - dni : float or Series + dni : numeric Direct Normal Irradiance Returns ------- - Series + beam : numeric + Beam component """ beam = dni * aoi_projection(surface_tilt, surface_azimuth, solar_zenith, solar_azimuth) - beam[beam < 0] = 0 + beam = np.maximum(beam, 0) return beam @@ -320,8 +323,7 @@ def total_irrad(surface_tilt, surface_azimuth, model='isotropic', model_perez='allsitescomposite1990', **kwargs): r""" - Determine diffuse irradiance from the sky on a - tilted surface. + Determine diffuse irradiance from the sky on a tilted surface. .. math:: @@ -331,25 +333,25 @@ def total_irrad(surface_tilt, surface_azimuth, Parameters ---------- - surface_tilt : float or Series. + surface_tilt : numeric Panel tilt from horizontal. - surface_azimuth : float or Series. + surface_azimuth : numeric Panel azimuth from north. - solar_zenith : float or Series. + solar_zenith : numeric Solar zenith angle. - solar_azimuth : float or Series. + solar_azimuth : numeric Solar azimuth angle. - dni : float or Series + dni : numeric Direct Normal Irradiance - ghi : float or Series + ghi : numeric Global horizontal irradiance - dhi : float or Series + dhi : numeric Diffuse horizontal irradiance - dni_extra : float or Series + dni_extra : numeric Extraterrestrial direct normal irradiance - airmass : float or Series + airmass : numeric Airmass - albedo : float + albedo : numeric Surface albedo surface_type : String Surface type. See grounddiffuse. @@ -360,8 +362,8 @@ def total_irrad(surface_tilt, surface_azimuth, Returns ------- - irradiance: DataFrame - Contains columns ``'poa_global', 'poa_direct', + irradiance : OrderedDict or DataFrame + Contains keys/columns ``'poa_global', 'poa_direct', 'poa_sky_diffuse', 'poa_ground_diffuse'``. """ @@ -399,47 +401,50 @@ def total_irrad(surface_tilt, surface_azimuth, diffuse = sky + ground total = beam + diffuse - all_irrad = pd.DataFrame() + all_irrad = OrderedDict() all_irrad['poa_global'] = total all_irrad['poa_direct'] = beam all_irrad['poa_diffuse'] = diffuse all_irrad['poa_sky_diffuse'] = sky all_irrad['poa_ground_diffuse'] = ground + if isinstance(total, pd.Series): + all_irrad = pd.DataFrame(all_irrad) + return all_irrad -# ToDo: keep this or not? wholmgren, 2014-11-03 def globalinplane(aoi, dni, poa_sky_diffuse, poa_ground_diffuse): r''' Determine the three components on in-plane irradiance - Combines in-plane irradaince compoents from the chosen diffuse translation, - ground reflection and beam irradiance algorithms into the total in-plane - irradiance. + Combines in-plane irradaince compoents from the chosen diffuse + translation, ground reflection and beam irradiance algorithms into + the total in-plane irradiance. Parameters ---------- + aoi : numeric + Angle of incidence of solar rays with respect to the module + surface, from :func:`aoi`. - aoi : float or Series - Angle of incidence of solar rays with respect - to the module surface, from :func:`aoi`. - - dni : float or Series - Direct normal irradiance (W/m^2), as measured - from a TMY file or calculated with a clearsky model. + dni : numeric + Direct normal irradiance (W/m^2), as measured from a TMY file or + calculated with a clearsky model. - poa_sky_diffuse : float or Series + poa_sky_diffuse : numeric Diffuse irradiance (W/m^2) in the plane of the modules, as calculated by a diffuse irradiance translation function - poa_ground_diffuse : float or Series + poa_ground_diffuse : numeric Ground reflected irradiance (W/m^2) in the plane of the modules, as calculated by an albedo model (eg. :func:`grounddiffuse`) Returns ------- - DataFrame with the following keys: + irrads : OrderedDict or DataFrame + Contains the following keys: + * ``poa_global`` : Total in-plane irradiance (W/m^2) * ``poa_direct`` : Total in-plane beam irradiance (W/m^2) * ``poa_diffuse`` : Total in-plane diffuse irradiance (W/m^2) @@ -450,13 +455,19 @@ def globalinplane(aoi, dni, poa_sky_diffuse, poa_ground_diffuse): :math:`< 0^{\circ}` is set to zero. ''' - poa_direct = pd.Series(dni * np.cos(np.radians(aoi))).clip_lower(0) + poa_direct = np.maximum(dni * np.cos(np.radians(aoi)), 0) poa_global = poa_direct + poa_sky_diffuse + poa_ground_diffuse poa_diffuse = poa_sky_diffuse + poa_ground_diffuse - return pd.DataFrame({'poa_global': poa_global, - 'poa_direct': poa_direct, - 'poa_diffuse': poa_diffuse}) + irrads = OrderedDict() + irrads['poa_global'] = poa_global + irrads['poa_direct'] = poa_direct + irrads['poa_diffuse'] = poa_diffuse + + if isinstance(poa_direct, pd.Series): + irrads = pd.DataFrame(irrads) + + return irrads def grounddiffuse(surface_tilt, ghi, albedo=.25, surface_type=None): @@ -464,42 +475,39 @@ def grounddiffuse(surface_tilt, ghi, albedo=.25, surface_type=None): Estimate diffuse irradiance from ground reflections given irradiance, albedo, and surface tilt - Function to determine the portion of irradiance on a tilted surface due - to ground reflections. Any of the inputs may be DataFrames or scalars. + Function to determine the portion of irradiance on a tilted surface + due to ground reflections. Any of the inputs may be DataFrames or + scalars. Parameters ---------- - surface_tilt : float or DataFrame - Surface tilt angles in decimal degrees. - SurfTilt must be >=0 and <=180. The tilt angle is defined as - degrees from horizontal (e.g. surface facing up = 0, surface facing - horizon = 90). + surface_tilt : numeric + Surface tilt angles in decimal degrees. Tilt must be >=0 and + <=180. The tilt angle is defined as degrees from horizontal + (e.g. surface facing up = 0, surface facing horizon = 90). - ghi : float or DataFrame + ghi : numeric Global horizontal irradiance in W/m^2. - albedo : float or DataFrame - Ground reflectance, typically 0.1-0.4 for - surfaces on Earth (land), may increase over snow, ice, etc. May also - be known as the reflection coefficient. Must be >=0 and <=1. - Will be overridden if surface_type is supplied. + albedo : numeric + Ground reflectance, typically 0.1-0.4 for surfaces on Earth + (land), may increase over snow, ice, etc. May also be known as + the reflection coefficient. Must be >=0 and <=1. Will be + overridden if surface_type is supplied. - surface_type: None or string in - ``'urban', 'grass', 'fresh grass', 'snow', 'fresh snow', - 'asphalt', 'concrete', 'aluminum', 'copper', - 'fresh steel', 'dirty steel'``. - Overrides albedo. + surface_type: None or string + If not None, overrides albedo. String can be one of ``'urban', + 'grass', 'fresh grass', 'snow', 'fresh snow', 'asphalt', + 'concrete', 'aluminum', 'copper', 'fresh steel', 'dirty steel'``. Returns ------- - - float or DataFrame - Ground reflected irradiances in W/m^2. + grounddiffuse : numeric + Ground reflected irradiances in W/m^2. References ---------- - [1] Loutzenhiser P.G. et. al. "Empirical validation of models to compute solar irradiance on inclined surfaces for building energy simulation" 2007, Solar Energy vol. 81. pp. 254-267. @@ -508,10 +516,10 @@ def grounddiffuse(surface_tilt, ghi, albedo=.25, surface_type=None): [2] albedos from: http://pvpmc.org/modeling-steps/incident-irradiance/plane-of-array-poa-irradiance/calculating-poa-irradiance/poa-ground-reflected/albedo/ - and http://en.wikipedia.org/wiki/Albedo ''' + pvl_logger.debug('diffuse_ground.get_diffuse_ground()') if surface_type is not None: @@ -531,55 +539,42 @@ def grounddiffuse(surface_tilt, ghi, albedo=.25, surface_type=None): def isotropic(surface_tilt, dhi): r''' - Determine diffuse irradiance from the sky on a - tilted surface using the isotropic sky model. + Determine diffuse irradiance from the sky on a tilted surface using + the isotropic sky model. .. math:: I_{d} = DHI \frac{1 + \cos\beta}{2} - Hottel and Woertz's model treats the sky as a uniform source of diffuse - irradiance. Thus the diffuse irradiance from the sky (ground reflected - irradiance is not included in this algorithm) on a tilted surface can - be found from the diffuse horizontal irradiance and the tilt angle of - the surface. + Hottel and Woertz's model treats the sky as a uniform source of + diffuse irradiance. Thus the diffuse irradiance from the sky (ground + reflected irradiance is not included in this algorithm) on a tilted + surface can be found from the diffuse horizontal irradiance and the + tilt angle of the surface. Parameters ---------- + surface_tilt : numeric + Surface tilt angle in decimal degrees. Tilt must be >=0 and + <=180. The tilt angle is defined as degrees from horizontal + (e.g. surface facing up = 0, surface facing horizon = 90) - surface_tilt : float or Series - Surface tilt angle in decimal degrees. - surface_tilt must be >=0 and <=180. The tilt angle is defined as - degrees from horizontal (e.g. surface facing up = 0, surface facing - horizon = 90) - - dhi : float or Series - Diffuse horizontal irradiance in W/m^2. - DHI must be >=0. - + dhi : numeric + Diffuse horizontal irradiance in W/m^2. DHI must be >=0. Returns ------- - float or Series - - The diffuse component of the solar radiation on an - arbitrarily tilted surface defined by the isotropic sky model as - given in Loutzenhiser et. al (2007) equation 3. - SkyDiffuse is the diffuse component ONLY and does not include the ground - reflected irradiance or the irradiance due to the beam. - SkyDiffuse is a column vector vector with a number of elements equal to - the input vector(s). - + diffuse : numeric + The sky diffuse component of the solar radiation. References ---------- + [1] Loutzenhiser P.G. et. al. "Empirical validation of models to + compute solar irradiance on inclined surfaces for building energy + simulation" 2007, Solar Energy vol. 81. pp. 254-267 - [1] Loutzenhiser P.G. et. al. "Empirical validation of models to compute - solar irradiance on inclined surfaces for building energy simulation" - 2007, Solar Energy vol. 81. pp. 254-267 - - [2] Hottel, H.C., Woertz, B.B., 1942. Evaluation of flat-plate solar heat - collector. Trans. ASME 64, 91. + [2] Hottel, H.C., Woertz, B.B., 1942. Evaluation of flat-plate solar + heat collector. Trans. ASME 64, 91. ''' pvl_logger.debug('diffuse_sky.isotropic()') @@ -589,8 +584,8 @@ def isotropic(surface_tilt, dhi): return sky_diffuse -def klucher(surface_tilt, surface_azimuth, dhi, ghi, - solar_zenith, solar_azimuth): +def klucher(surface_tilt, surface_azimuth, dhi, ghi, solar_zenith, + solar_azimuth): r''' Determine diffuse irradiance from the sky on a tilted surface using Klucher's 1979 model @@ -610,53 +605,40 @@ def klucher(surface_tilt, surface_azimuth, dhi, ghi, (ground reflected irradiance is not included in this algorithm) on a tilted surface using the surface tilt angle, surface azimuth angle, diffuse horizontal irradiance, direct normal irradiance, global - horizontal irradiance, extraterrestrial irradiance, sun zenith angle, - and sun azimuth angle. + horizontal irradiance, extraterrestrial irradiance, sun zenith + angle, and sun azimuth angle. Parameters ---------- + surface_tilt : numeric + Surface tilt angles in decimal degrees. surface_tilt must be >=0 + and <=180. The tilt angle is defined as degrees from horizontal + (e.g. surface facing up = 0, surface facing horizon = 90) - surface_tilt : float or Series - Surface tilt angles in decimal degrees. - surface_tilt must be >=0 and <=180. The tilt angle is defined as - degrees from horizontal (e.g. surface facing up = 0, surface facing - horizon = 90) - - surface_azimuth : float or Series - Surface azimuth angles in decimal degrees. - surface_azimuth must be >=0 and <=360. The Azimuth convention is defined - as degrees east of north (e.g. North = 0, South=180 East = 90, - West = 270). + surface_azimuth : numeric + Surface azimuth angles in decimal degrees. surface_azimuth must + be >=0 and <=360. The Azimuth convention is defined as degrees + east of north (e.g. North = 0, South=180 East = 90, West = 270). - dhi : float or Series - diffuse horizontal irradiance in W/m^2. - DHI must be >=0. + dhi : numeric + Diffuse horizontal irradiance in W/m^2. DHI must be >=0. - ghi : float or Series - Global irradiance in W/m^2. - DNI must be >=0. + ghi : numeric + Global irradiance in W/m^2. DNI must be >=0. - solar_zenith : float or Series - apparent (refraction-corrected) zenith - angles in decimal degrees. - solar_zenith must be >=0 and <=180. + solar_zenith : numeric + Apparent (refraction-corrected) zenith angles in decimal + degrees. solar_zenith must be >=0 and <=180. - solar_azimuth : float or Series - Sun azimuth angles in decimal degrees. - solar_azimuth must be >=0 and <=360. The Azimuth convention is defined - as degrees east of north (e.g. North = 0, East = 90, West = 270). + solar_azimuth : numeric + Sun azimuth angles in decimal degrees. solar_azimuth must be >=0 + and <=360. The Azimuth convention is defined as degrees east of + north (e.g. North = 0, East = 90, West = 270). Returns ------- - float or Series. - - The diffuse component of the solar radiation on an - arbitrarily tilted surface defined by the Klucher model as given in - Loutzenhiser et. al (2007) equation 4. - SkyDiffuse is the diffuse component ONLY and does not include the ground - reflected irradiance or the irradiance due to the beam. - SkyDiffuse is a column vector vector with a number of elements equal to - the input vector(s). + diffuse : numeric + The sky diffuse component of the solar radiation. References ---------- @@ -693,78 +675,70 @@ def klucher(surface_tilt, surface_azimuth, dhi, ghi, def haydavies(surface_tilt, surface_azimuth, dhi, dni, dni_extra, solar_zenith=None, solar_azimuth=None, projection_ratio=None): r''' - Determine diffuse irradiance from the sky on a - tilted surface using Hay & Davies' 1980 model + Determine diffuse irradiance from the sky on a tilted surface using + Hay & Davies' 1980 model .. math:: I_{d} = DHI ( A R_b + (1 - A) (\frac{1 + \cos\beta}{2}) ) - Hay and Davies' 1980 model determines the diffuse irradiance from the sky - (ground reflected irradiance is not included in this algorithm) on a - tilted surface using the surface tilt angle, surface azimuth angle, - diffuse horizontal irradiance, direct normal irradiance, - extraterrestrial irradiance, sun zenith angle, and sun azimuth angle. + Hay and Davies' 1980 model determines the diffuse irradiance from + the sky (ground reflected irradiance is not included in this + algorithm) on a tilted surface using the surface tilt angle, surface + azimuth angle, diffuse horizontal irradiance, direct normal + irradiance, extraterrestrial irradiance, sun zenith angle, and sun + azimuth angle. Parameters ---------- + surface_tilt : numeric + Surface tilt angles in decimal degrees. The tilt angle is + defined as degrees from horizontal (e.g. surface facing up = 0, + surface facing horizon = 90) - surface_tilt : float or Series - Surface tilt angles in decimal degrees. - The tilt angle is defined as - degrees from horizontal (e.g. surface facing up = 0, surface facing - horizon = 90) - - surface_azimuth : float or Series - Surface azimuth angles in decimal degrees. - The azimuth convention is defined - as degrees east of north (e.g. North=0, South=180, East=90, - West=270). + surface_azimuth : numeric + Surface azimuth angles in decimal degrees. The azimuth + convention is defined as degrees east of north (e.g. North=0, + South=180, East=90, West=270). - dhi : float or Series + dhi : numeric Diffuse horizontal irradiance in W/m^2. - dni : float or Series + dni : numeric Direct normal irradiance in W/m^2. - dni_extra : float or Series + dni_extra : numeric Extraterrestrial normal irradiance in W/m^2. - solar_zenith : None, float or Series - Solar apparent (refraction-corrected) zenith - angles in decimal degrees. - Must supply ``solar_zenith`` and ``solar_azimuth`` or supply - ``projection_ratio``. + solar_zenith : None or numeric + Solar apparent (refraction-corrected) zenith angles in decimal + degrees. Must supply ``solar_zenith`` and ``solar_azimuth`` or + supply ``projection_ratio``. - solar_azimuth : None, float or Series - Solar azimuth angles in decimal degrees. - Must supply ``solar_zenith`` and ``solar_azimuth`` or supply + solar_azimuth : None or numeric + Solar azimuth angles in decimal degrees. Must supply + ``solar_zenith`` and ``solar_azimuth`` or supply ``projection_ratio``. - projection_ratio : None, float or Series - Ratio of angle of incidence projection to solar zenith - angle projection. - Must supply ``solar_zenith`` and ``solar_azimuth`` or supply - ``projection_ratio``. + projection_ratio : None or numeric + Ratio of angle of incidence projection to solar zenith angle + projection. Must supply ``solar_zenith`` and ``solar_azimuth`` + or supply ``projection_ratio``. Returns -------- - - sky_diffuse : float or Series - The diffuse component of the solar radiation on an - arbitrarily tilted surface defined by the Perez model as given in - reference [3]. Does not include the - ground reflected irradiance or the irradiance due to the beam. + sky_diffuse : numeric + The sky diffuse component of the solar radiation. References ----------- - [1] Loutzenhiser P.G. et. al. "Empirical validation of models to compute - solar irradiance on inclined surfaces for building energy simulation" - 2007, Solar Energy vol. 81. pp. 254-267 + [1] Loutzenhiser P.G. et. al. "Empirical validation of models to + compute solar irradiance on inclined surfaces for building energy + simulation" 2007, Solar Energy vol. 81. pp. 254-267 - [2] Hay, J.E., Davies, J.A., 1980. Calculations of the solar radiation - incident on an inclined surface. In: Hay, J.E., Won, T.K. (Eds.), Proc. of - First Canadian Solar Radiation Data Workshop, 59. Ministry of Supply - and Services, Canada. + [2] Hay, J.E., Davies, J.A., 1980. Calculations of the solar + radiation incident on an inclined surface. In: Hay, J.E., Won, T.K. + (Eds.), Proc. of First Canadian Solar Radiation Data Workshop, 59. + Ministry of Supply and Services, Canada. ''' pvl_logger.debug('diffuse_sky.haydavies()') @@ -781,12 +755,12 @@ def haydavies(surface_tilt, surface_azimuth, dhi, dni, dni_extra, # Anisotropy Index AI = dni / dni_extra - # these are actually the () and [] sub-terms of the second term of eqn 7 + # these are the () and [] sub-terms of the second term of eqn 7 term1 = 1 - AI term2 = 0.5 * (1 + tools.cosd(surface_tilt)) sky_diffuse = dhi * (AI * Rb + term1 * term2) - sky_diffuse[sky_diffuse < 0] = 0 + sky_diffuse = np.maximum(sky_diffuse, 0) return sky_diffuse @@ -794,8 +768,8 @@ def haydavies(surface_tilt, surface_azimuth, dhi, dni, dni_extra, def reindl(surface_tilt, surface_azimuth, dhi, dni, ghi, dni_extra, solar_zenith, solar_azimuth): r''' - Determine diffuse irradiance from the sky on a - tilted surface using Reindl's 1990 model + Determine diffuse irradiance from the sky on a tilted surface using + Reindl's 1990 model .. math:: @@ -806,62 +780,48 @@ def reindl(surface_tilt, surface_azimuth, dhi, dni, ghi, dni_extra, (ground reflected irradiance is not included in this algorithm) on a tilted surface using the surface tilt angle, surface azimuth angle, diffuse horizontal irradiance, direct normal irradiance, global - horizontal irradiance, extraterrestrial irradiance, sun zenith angle, - and sun azimuth angle. + horizontal irradiance, extraterrestrial irradiance, sun zenith + angle, and sun azimuth angle. Parameters ---------- + surface_tilt : numeric + Surface tilt angles in decimal degrees. The tilt angle is + defined as degrees from horizontal (e.g. surface facing up = 0, + surface facing horizon = 90) - surface_tilt : float or Series. - Surface tilt angles in decimal degrees. - The tilt angle is defined as - degrees from horizontal (e.g. surface facing up = 0, surface facing - horizon = 90) - - surface_azimuth : float or Series. - Surface azimuth angles in decimal degrees. - The Azimuth convention is defined - as degrees east of north (e.g. North = 0, South=180 East = 90, - West = 270). + surface_azimuth : numeric + Surface azimuth angles in decimal degrees. The azimuth + convention is defined as degrees east of north (e.g. North = 0, + South=180 East = 90, West = 270). - dhi : float or Series. + dhi : numeric diffuse horizontal irradiance in W/m^2. - dni : float or Series. + dni : numeric direct normal irradiance in W/m^2. - ghi: float or Series. + ghi: numeric Global irradiance in W/m^2. - dni_extra : float or Series. - extraterrestrial normal irradiance in W/m^2. + dni_extra : numeric + Extraterrestrial normal irradiance in W/m^2. - solar_zenith : float or Series. - apparent (refraction-corrected) zenith - angles in decimal degrees. + solar_zenith : numeric + Apparent (refraction-corrected) zenith angles in decimal degrees. - solar_azimuth : float or Series. - Sun azimuth angles in decimal degrees. - The Azimuth convention is defined - as degrees east of north (e.g. North = 0, East = 90, West = 270). + solar_azimuth : numeric + Sun azimuth angles in decimal degrees. The azimuth convention is + defined as degrees east of north (e.g. North = 0, East = 90, + West = 270). Returns ------- - - poa_sky_diffuse : float or Series. - - The diffuse component of the solar radiation on an - arbitrarily tilted surface defined by the Reindl model as given in - Loutzenhiser et. al (2007) equation 8. - SkyDiffuse is the diffuse component ONLY and does not include the - ground reflected irradiance or the irradiance due to the beam. - SkyDiffuse is a column vector vector with a number of elements equal to - the input vector(s). - + poa_sky_diffuse : numeric + The sky diffuse component of the solar radiation. Notes ----- - The poa_sky_diffuse calculation is generated from the Loutzenhiser et al. (2007) paper, equation 8. Note that I have removed the beam and ground reflectance portion of the equation and this generates ONLY the diffuse @@ -870,16 +830,15 @@ def reindl(surface_tilt, surface_azimuth, dhi, dni, ghi, dni_extra, References ---------- + [1] Loutzenhiser P.G. et. al. "Empirical validation of models to + compute solar irradiance on inclined surfaces for building energy + simulation" 2007, Solar Energy vol. 81. pp. 254-267 - [1] Loutzenhiser P.G. et. al. "Empirical validation of models to compute - solar irradiance on inclined surfaces for building energy simulation" - 2007, Solar Energy vol. 81. pp. 254-267 - - [2] Reindl, D.T., Beckmann, W.A., Duffie, J.A., 1990a. Diffuse fraction - correlations. Solar Energy 45(1), 1-7. + [2] Reindl, D.T., Beckmann, W.A., Duffie, J.A., 1990a. Diffuse + fraction correlations. Solar Energy 45(1), 1-7. - [3] Reindl, D.T., Beckmann, W.A., Duffie, J.A., 1990b. Evaluation of hourly - tilted surface radiation models. Solar Energy 45(1), 9-17. + [3] Reindl, D.T., Beckmann, W.A., Duffie, J.A., 1990b. Evaluation of + hourly tilted surface radiation models. Solar Energy 45(1), 9-17. ''' pvl_logger.debug('diffuse_sky.reindl()') @@ -897,56 +856,51 @@ def reindl(surface_tilt, surface_azimuth, dhi, dni, ghi, dni_extra, # DNI projected onto horizontal HB = dni * cos_solar_zenith - HB[HB < 0] = 0 + HB = np.maximum(HB, 0) - # these are actually the () and [] sub-terms of the second term of eqn 8 + # these are the () and [] sub-terms of the second term of eqn 8 term1 = 1 - AI term2 = 0.5 * (1 + tools.cosd(surface_tilt)) term3 = 1 + np.sqrt(HB / ghi) * (tools.sind(0.5 * surface_tilt) ** 3) sky_diffuse = dhi * (AI * Rb + term1 * term2 * term3) - sky_diffuse[sky_diffuse < 0] = 0 + sky_diffuse = np.maximum(sky_diffuse, 0) return sky_diffuse def king(surface_tilt, dhi, ghi, solar_zenith): ''' - Determine diffuse irradiance from the sky on a tilted surface using the - King model. + Determine diffuse irradiance from the sky on a tilted surface using + the King model. - King's model determines the diffuse irradiance from the sky - (ground reflected irradiance is not included in this algorithm) on a - tilted surface using the surface tilt angle, diffuse horizontal - irradiance, global horizontal irradiance, and sun zenith angle. Note - that this model is not well documented and has not been published in - any fashion (as of January 2012). + King's model determines the diffuse irradiance from the sky (ground + reflected irradiance is not included in this algorithm) on a tilted + surface using the surface tilt angle, diffuse horizontal irradiance, + global horizontal irradiance, and sun zenith angle. Note that this + model is not well documented and has not been published in any + fashion (as of January 2012). Parameters ---------- + surface_tilt : numeric + Surface tilt angles in decimal degrees. The tilt angle is + defined as degrees from horizontal (e.g. surface facing up = 0, + surface facing horizon = 90) - surface_tilt : float or Series - Surface tilt angles in decimal degrees. - The tilt angle is defined as - degrees from horizontal (e.g. surface facing up = 0, surface facing - horizon = 90) - - dhi : float or Series + dhi : numeric Diffuse horizontal irradiance in W/m^2. - ghi : float or Series + ghi : numeric Global horizontal irradiance in W/m^2. - solar_zenith : float or Series - Apparent (refraction-corrected) zenith - angles in decimal degrees. + solar_zenith : numeric + Apparent (refraction-corrected) zenith angles in decimal degrees. Returns -------- - poa_sky_diffuse : float or Series - The diffuse component of the solar radiation on an - arbitrarily tilted surface as given by a model developed by - David L. King at Sandia National Laboratories. + poa_sky_diffuse : numeric + The diffuse component of the solar radiation. ''' pvl_logger.debug('diffuse_sky.king()') @@ -954,7 +908,7 @@ def king(surface_tilt, dhi, ghi, solar_zenith): sky_diffuse = (dhi * ((1 + tools.cosd(surface_tilt))) / 2 + ghi * ((0.012 * solar_zenith - 0.04)) * ((1 - tools.cosd(surface_tilt))) / 2) - sky_diffuse[sky_diffuse < 0] = 0 + sky_diffuse = np.maximum(sky_diffuse, 0) return sky_diffuse @@ -974,42 +928,37 @@ def perez(surface_tilt, surface_azimuth, dhi, dni, dni_extra, pressure-corrected) airmass. Optionally a selector may be used to use any of Perez's model coefficient sets. - Parameters ---------- - - surface_tilt : float or Series + surface_tilt : numeric Surface tilt angles in decimal degrees. surface_tilt must be >=0 and <=180. The tilt angle is defined as degrees from horizontal (e.g. surface facing up = 0, surface facing horizon = 90) - surface_azimuth : float or Series + surface_azimuth : numeric Surface azimuth angles in decimal degrees. surface_azimuth must - be >=0 and <=360. The Azimuth convention is defined as degrees + be >=0 and <=360. The azimuth convention is defined as degrees east of north (e.g. North = 0, South=180 East = 90, West = 270). - dhi : float or Series - Diffuse horizontal irradiance in W/m^2. - DHI must be >=0. + dhi : numeric + Diffuse horizontal irradiance in W/m^2. DHI must be >=0. - dni : float or Series - Direct normal irradiance in W/m^2. - DNI must be >=0. + dni : numeric + Direct normal irradiance in W/m^2. DNI must be >=0. - dni_extra : float or Series + dni_extra : numeric Extraterrestrial normal irradiance in W/m^2. - solar_zenith : float or Series - apparent (refraction-corrected) zenith - angles in decimal degrees. - solar_zenith must be >=0 and <=180. + solar_zenith : numeric + apparent (refraction-corrected) zenith angles in decimal + degrees. solar_zenith must be >=0 and <=180. - solar_azimuth : float or Series + solar_azimuth : numeric Sun azimuth angles in decimal degrees. solar_azimuth must be >=0 - and <=360. The Azimuth convention is defined as degrees east of + and <=360. The azimuth convention is defined as degrees east of north (e.g. North = 0, East = 90, West = 270). - airmass : float or Series + airmass : numeric Relative (not pressure-corrected) airmass values. If AM is a DataFrame it must be of the same size as all other DataFrame inputs. AM must be >=0 (careful using the 1/sec(z) model of AM @@ -1035,19 +984,12 @@ def perez(surface_tilt, surface_azimuth, dhi, dni, dni_extra, Returns -------- - - float or Series - - The diffuse component of the solar radiation on an arbitrarily - tilted surface defined by the Perez model as given in reference - [3]. SkyDiffuse is the diffuse component ONLY and does not - include the ground reflected irradiance or the irradiance due to - the beam. - + sky_diffuse : numeric + The sky diffuse component of the solar radiation on a tilted + surface. Array input is currently converted to Series output. References ---------- - [1] Loutzenhiser P.G. et. al. "Empirical validation of models to compute solar irradiance on inclined surfaces for building energy simulation" 2007, Solar Energy vol. 81. pp. 254-267 @@ -1067,15 +1009,22 @@ def perez(surface_tilt, surface_azimuth, dhi, dni, dni_extra, kappa = 1.041 # for solar_zenith in radians z = np.radians(solar_zenith) # convert to radians + # delta is the sky's "brightness" + delta = dhi * airmass / dni_extra + # epsilon is the sky's "clearness" eps = ((dhi + dni) / dhi + kappa * (z ** 3)) / (1 + kappa * (z ** 3)) + # numpy indexing below will not work with a Series + if isinstance(eps, pd.Series): + eps = eps.values + # Perez et al define clearness bins according to the following # rules. 1 = overcast ... 8 = clear (these names really only make # sense for small zenith angles, but...) these values will # eventually be used as indicies for coeffecient look ups - ebin = eps.copy() - ebin[(eps < 1.065)] = 1 + ebin = np.zeros_like(eps, dtype=np.int8) + ebin[eps < 1.065] = 1 ebin[(eps >= 1.065) & (eps < 1.23)] = 2 ebin[(eps >= 1.23) & (eps < 1.5)] = 3 ebin[(eps >= 1.5) & (eps < 1.95)] = 4 @@ -1084,144 +1033,477 @@ def perez(surface_tilt, surface_azimuth, dhi, dni, dni_extra, ebin[(eps >= 4.5) & (eps < 6.2)] = 7 ebin[eps >= 6.2] = 8 - ebin = ebin - 1 # correct for 0 indexing in coeffecient lookup - - # remove night time values - ebin = ebin.dropna().astype(int) - - # This is added because in cases where the sun is below the horizon - # (var.solar_zenith > 90) but there is still diffuse horizontal - # light (var.DHI>0), it is possible that the airmass (var.AM) could - # be NaN, which messes up later calculations. Instead, if the sun is - # down, and there is still var.DHI, we set the airmass to the - # airmass value on the horizon (approximately 37-38). - # var.AM(var.solar_zenith >=90 & var.DHI >0) = 37; - - # var.DNI_ET[var.DNI_ET==0] = .00000001 #very hacky, fix this - - # delta is the sky's "brightness" - delta = dhi * airmass / dni_extra + # correct for 0 indexing in coeffecient lookup + # later, ebin = -1 will yield nan coefficients + ebin -= 1 # The various possible sets of Perez coefficients are contained # in a subfunction to clean up the code. F1c, F2c = _get_perez_coefficients(model) - F1 = (F1c[ebin, 0] + F1c[ebin, 1] * delta[ebin.index] + - F1c[ebin, 2] * z[ebin.index]) - F1[F1 < 0] = 0 - F1 = F1.astype(float) + # results in invalid eps (ebin = -1) being mapped to nans + nans = np.array([np.nan, np.nan, np.nan]) + F1c = np.vstack((F1c, nans)) + F2c = np.vstack((F2c, nans)) - F2 = (F2c[ebin, 0] + F2c[ebin, 1] * delta[ebin.index] + - F2c[ebin, 2] * z[ebin.index]) - F2[F2 < 0] = 0 - F2 = F2.astype(float) + F1 = (F1c[ebin, 0] + F1c[ebin, 1] * delta + F1c[ebin, 2] * z) + F1 = np.maximum(F1, 0) + + F2 = (F2c[ebin, 0] + F2c[ebin, 1] * delta + F2c[ebin, 2] * z) + F2 = np.maximum(F2, 0) A = aoi_projection(surface_tilt, surface_azimuth, solar_zenith, solar_azimuth) - A[A < 0] = 0 + A = np.maximum(A, 0) B = tools.cosd(solar_zenith) - B[B < tools.cosd(85)] = tools.cosd(85) + B = np.maximum(B, tools.cosd(85)) # Calculate Diffuse POA from sky dome - term1 = 0.5 * (1 - F1) * (1 + tools.cosd(surface_tilt)) - term2 = F1 * A[ebin.index] / B[ebin.index] + term2 = F1 * A / B term3 = F2 * tools.sind(surface_tilt) - sky_diffuse = dhi * (term1 + term2 + term3) - sky_diffuse[sky_diffuse < 0] = 0 - sky_diffuse[airmass.isnull()] = 0 + sky_diffuse = np.maximum(dhi * (term1 + term2 + term3), 0) + + # we've preserved the input type until now, so don't ruin it! + if isinstance(sky_diffuse, pd.Series): + sky_diffuse[np.isnan(airmass)] = 0 + else: + sky_diffuse = np.where(np.isnan(airmass), 0, sky_diffuse) return sky_diffuse -def _get_perez_coefficients(perezmodel): - ''' - Find coefficients for the Perez model +def disc(ghi, zenith, times, pressure=101325): + """ + Estimate Direct Normal Irradiance from Global Horizontal Irradiance + using the DISC model. + + The DISC algorithm converts global horizontal irradiance to direct + normal irradiance through empirical relationships between the global + and direct clearness indices. Parameters ---------- + ghi : numeric + Global horizontal irradiance in W/m^2. - perezmodel : string (optional, default='allsitescomposite1990') - - a character string which selects the desired set of Perez - coefficients. If model is not provided as an input, the default, - '1990' will be used. + solar_zenith : numeric + True (not refraction-corrected) solar zenith angles in decimal + degrees. - All possible model selections are: + times : DatetimeIndex - * '1990' - * 'allsitescomposite1990' (same as '1990') - * 'allsitescomposite1988' - * 'sandiacomposite1988' - * 'usacomposite1988' - * 'france1988' - * 'phoenix1988' - * 'elmonte1988' - * 'osage1988' - * 'albuquerque1988' - * 'capecanaveral1988' - * 'albany1988' + pressure : numeric + Site pressure in Pascal. Returns + ------- + output : OrderedDict or DataFrame + Contains the following keys: + + * ``dni``: The modeled direct normal irradiance + in W/m^2 provided by the + Direct Insolation Simulation Code (DISC) model. + * ``kt``: Ratio of global to extraterrestrial + irradiance on a horizontal plane. + * ``airmass``: Airmass + + References + ---------- + [1] Maxwell, E. L., "A Quasi-Physical Model for Converting Hourly + Global Horizontal to Direct Normal Insolation", Technical + Report No. SERI/TR-215-3087, Golden, CO: Solar Energy Research + Institute, 1987. + + [2] J.W. "Fourier series representation of the position of the sun". + Found at: + http://www.mail-archive.com/sundial@uni-koeln.de/msg01050.html on + January 12, 2012 + + See Also -------- + atmosphere.alt2pres + dirint + """ + + # in principle, the dni_extra calculation could be done by + # pvlib's function. However, this is the algorithm used in + # the DISC paper + + doy = times.dayofyear + + dayangle = 2. * np.pi*(doy - 1) / 365 + + re = (1.00011 + 0.034221*np.cos(dayangle) + 0.00128*np.sin(dayangle) + + 0.000719*np.cos(2.*dayangle) + 7.7e-5*np.sin(2.*dayangle)) + + I0 = re * 1370. + I0h = I0 * np.cos(np.radians(zenith)) + + am = atmosphere.relativeairmass(zenith, model='kasten1966') + am = atmosphere.absoluteairmass(am, pressure) + + kt = ghi / I0h + kt = np.maximum(kt, 0) + # powers of kt will be used repeatedly, so compute only once + kt2 = kt * kt # about the same as kt ** 2 + kt3 = kt2 * kt # 5-10x faster than kt ** 3 + + bools = (kt <= 0.6) + a = np.where(bools, + 0.512 - 1.56*kt + 2.286*kt2 - 2.222*kt3, + -5.743 + 21.77*kt - 27.49*kt2 + 11.56*kt3) + b = np.where(bools, + 0.37 + 0.962*kt, + 41.4 - 118.5*kt + 66.05*kt2 + 31.9*kt3) + c = np.where(bools, + -0.28 + 0.932*kt - 2.048*kt2, + -47.01 + 184.2*kt - 222.0*kt2 + 73.81*kt3) + + delta_kn = a + b * np.exp(c*am) + + Knc = 0.866 - 0.122*am + 0.0121*am**2 - 0.000653*am**3 + 1.4e-05*am**4 + Kn = Knc - delta_kn + + dni = Kn * I0 + + dni = np.where((zenith > 87) | (ghi < 0) | (dni < 0), 0, dni) + + output = OrderedDict() + output['dni'] = dni + output['kt'] = kt + output['airmass'] = am + + if isinstance(times, pd.DatetimeIndex): + output = pd.DataFrame(output, index=times) + + return output + + +def dirint(ghi, zenith, times, pressure=101325., use_delta_kt_prime=True, + temp_dew=None): + """ + Determine DNI from GHI using the DIRINT modification of the DISC + model. + + Implements the modified DISC model known as "DIRINT" introduced in + [1]. DIRINT predicts direct normal irradiance (DNI) from measured + global horizontal irradiance (GHI). DIRINT improves upon the DISC + model by using time-series GHI data and dew point temperature + information. The effectiveness of the DIRINT model improves with + each piece of information provided. + + Parameters + ---------- + ghi : array-like + Global horizontal irradiance in W/m^2. + + zenith : array-like + True (not refraction-corrected) zenith angles in decimal + degrees. If Z is a vector it must be of the same size as all + other vector inputs. Z must be >=0 and <=180. + + times : DatetimeIndex + + pressure : float or array-like + The site pressure in Pascal. Pressure may be measured or an + average pressure may be calculated from site altitude. + + use_delta_kt_prime : bool + Indicates if the user would like to utilize the time-series + nature of the GHI measurements. A value of ``False`` will not + use the time-series improvements, any other numeric value will + use time-series improvements. It is recommended that time-series + data only be used if the time between measured data points is + less than 1.5 hours. If none of the input arguments are vectors, + then time-series improvements are not used (because it's not a + time-series). If True, input data must be Series. + + temp_dew : None, float, or array-like + Surface dew point temperatures, in degrees C. Values of temp_dew + may be numeric or NaN. Any single time period point with a + DewPtTemp=NaN does not have dew point improvements applied. If + DewPtTemp is not provided, then dew point improvements are not + applied. + + Returns + ------- + dni : array-like + The modeled direct normal irradiance in W/m^2 provided by the + DIRINT model. - F1coeffs : array - F1 coefficients for the Perez model - F2coeffs : array - F2 coefficients for the Perez model + Notes + ----- + DIRINT model requires time series data (ie. one of the inputs must + be a vector of length > 2. References ---------- + [1] Perez, R., P. Ineichen, E. Maxwell, R. Seals and A. Zelenka, + (1992). "Dynamic Global-to-Direct Irradiance Conversion Models". + ASHRAE Transactions-Research Series, pp. 354-369 - [1] Loutzenhiser P.G. et. al. "Empirical validation of models to compute - solar irradiance on inclined surfaces for building energy simulation" - 2007, Solar Energy vol. 81. pp. 254-267 + [2] Maxwell, E. L., "A Quasi-Physical Model for Converting Hourly + Global Horizontal to Direct Normal Insolation", Technical Report No. + SERI/TR-215-3087, Golden, CO: Solar Energy Research Institute, 1987. + """ - [2] Perez, R., Seals, R., Ineichen, P., Stewart, R., Menicucci, D., 1987. - A new simplified version of the Perez diffuse irradiance model for tilted - surfaces. Solar Energy 39(3), 221-232. + disc_out = disc(ghi, zenith, times, pressure=pressure) + dni = disc_out['dni'] + kt = disc_out['kt'] + am = disc_out['airmass'] - [3] Perez, R., Ineichen, P., Seals, R., Michalsky, J., Stewart, R., 1990. - Modeling daylight availability and irradiance components from direct - and global irradiance. Solar Energy 44 (5), 271-289. + kt_prime = kt / (1.031 * np.exp(-1.4 / (0.9 + 9.4 / am)) + 0.1) + kt_prime = np.minimum(kt_prime, 0.82) # From SRRL code - [4] Perez, R. et. al 1988. "The Development and Verification of the - Perez Diffuse Radiation Model". SAND88-7030 + # wholmgren: + # the use_delta_kt_prime statement is a port of the MATLAB code. + # I am confused by the abs() in the delta_kt_prime calculation. + # It is not the absolute value of the central difference. + # current implementation requires that kt_prime is a Series + if use_delta_kt_prime: + delta_kt_prime = 0.5*((kt_prime - kt_prime.shift(1)).abs().add( + (kt_prime - kt_prime.shift(-1)).abs(), + fill_value=0)) + else: + delta_kt_prime = pd.Series(-1, index=times) - ''' - coeffdict = { - 'allsitescomposite1990': [ - [-0.0080, 0.5880, -0.0620, -0.0600, 0.0720, -0.0220], - [0.1300, 0.6830, -0.1510, -0.0190, 0.0660, -0.0290], - [0.3300, 0.4870, -0.2210, 0.0550, -0.0640, -0.0260], - [0.5680, 0.1870, -0.2950, 0.1090, -0.1520, -0.0140], - [0.8730, -0.3920, -0.3620, 0.2260, -0.4620, 0.0010], - [1.1320, -1.2370, -0.4120, 0.2880, -0.8230, 0.0560], - [1.0600, -1.6000, -0.3590, 0.2640, -1.1270, 0.1310], - [0.6780, -0.3270, -0.2500, 0.1560, -1.3770, 0.2510]], - 'allsitescomposite1988': [ - [-0.0180, 0.7050, -0.071, -0.0580, 0.1020, -0.0260], - [0.1910, 0.6450, -0.1710, 0.0120, 0.0090, -0.0270], - [0.4400, 0.3780, -0.2560, 0.0870, -0.1040, -0.0250], - [0.7560, -0.1210, -0.3460, 0.1790, -0.3210, -0.0080], - [0.9960, -0.6450, -0.4050, 0.2600, -0.5900, 0.0170], - [1.0980, -1.2900, -0.3930, 0.2690, -0.8320, 0.0750], - [0.9730, -1.1350, -0.3780, 0.1240, -0.2580, 0.1490], - [0.6890, -0.4120, -0.2730, 0.1990, -1.6750, 0.2370]], - 'sandiacomposite1988': [ - [-0.1960, 1.0840, -0.0060, -0.1140, 0.1800, -0.0190], - [0.2360, 0.5190, -0.1800, -0.0110, 0.0200, -0.0380], - [0.4540, 0.3210, -0.2550, 0.0720, -0.0980, -0.0460], - [0.8660, -0.3810, -0.3750, 0.2030, -0.4030, -0.0490], - [1.0260, -0.7110, -0.4260, 0.2730, -0.6020, -0.0610], - [0.9780, -0.9860, -0.3500, 0.2800, -0.9150, -0.0240], - [0.7480, -0.9130, -0.2360, 0.1730, -1.0450, 0.0650], - [0.3180, -0.7570, 0.1030, 0.0620, -1.6980, 0.2360]], - 'usacomposite1988': [ - [-0.0340, 0.6710, -0.0590, -0.0590, 0.0860, -0.0280], - [0.2550, 0.4740, -0.1910, 0.0180, -0.0140, -0.0330], + if temp_dew is not None: + w = pd.Series(np.exp(0.07 * temp_dew - 0.075), index=times) + else: + w = pd.Series(-1, index=times) + + # @wholmgren: the following bin assignments use MATLAB's 1-indexing. + # Later, we'll subtract 1 to conform to Python's 0-indexing. + + # Create kt_prime bins + kt_prime_bin = pd.Series(0, index=times, dtype=np.int64) + kt_prime_bin[(kt_prime >= 0) & (kt_prime < 0.24)] = 1 + kt_prime_bin[(kt_prime >= 0.24) & (kt_prime < 0.4)] = 2 + kt_prime_bin[(kt_prime >= 0.4) & (kt_prime < 0.56)] = 3 + kt_prime_bin[(kt_prime >= 0.56) & (kt_prime < 0.7)] = 4 + kt_prime_bin[(kt_prime >= 0.7) & (kt_prime < 0.8)] = 5 + kt_prime_bin[(kt_prime >= 0.8) & (kt_prime <= 1)] = 6 + + # Create zenith angle bins + zenith_bin = pd.Series(0, index=times, dtype=np.int64) + zenith_bin[(zenith >= 0) & (zenith < 25)] = 1 + zenith_bin[(zenith >= 25) & (zenith < 40)] = 2 + zenith_bin[(zenith >= 40) & (zenith < 55)] = 3 + zenith_bin[(zenith >= 55) & (zenith < 70)] = 4 + zenith_bin[(zenith >= 70) & (zenith < 80)] = 5 + zenith_bin[(zenith >= 80)] = 6 + + # Create the bins for w based on dew point temperature + w_bin = pd.Series(0, index=times, dtype=np.int64) + w_bin[(w >= 0) & (w < 1)] = 1 + w_bin[(w >= 1) & (w < 2)] = 2 + w_bin[(w >= 2) & (w < 3)] = 3 + w_bin[(w >= 3)] = 4 + w_bin[(w == -1)] = 5 + + # Create delta_kt_prime binning. + delta_kt_prime_bin = pd.Series(0, index=times, dtype=np.int64) + delta_kt_prime_bin[(delta_kt_prime >= 0) & (delta_kt_prime < 0.015)] = 1 + delta_kt_prime_bin[(delta_kt_prime >= 0.015) & + (delta_kt_prime < 0.035)] = 2 + delta_kt_prime_bin[(delta_kt_prime >= 0.035) & (delta_kt_prime < 0.07)] = 3 + delta_kt_prime_bin[(delta_kt_prime >= 0.07) & (delta_kt_prime < 0.15)] = 4 + delta_kt_prime_bin[(delta_kt_prime >= 0.15) & (delta_kt_prime < 0.3)] = 5 + delta_kt_prime_bin[(delta_kt_prime >= 0.3) & (delta_kt_prime <= 1)] = 6 + delta_kt_prime_bin[delta_kt_prime == -1] = 7 + + # get the coefficients + coeffs = _get_dirint_coeffs() + + # subtract 1 to account for difference between MATLAB-style bin + # assignment and Python-style array lookup. + dirint_coeffs = coeffs[kt_prime_bin-1, zenith_bin-1, + delta_kt_prime_bin-1, w_bin-1] + + # convert unassigned bins to nan + dirint_coeffs = np.where((kt_prime_bin == 0) | (zenith_bin == 0) | + (w_bin == 0) | (delta_kt_prime_bin == 0), + np.nan, dirint_coeffs) + + dni *= dirint_coeffs + + return dni + + +def erbs(ghi, zenith, doy): + r""" + Estimate DNI and DHI from GHI using the Erbs model. + + The Erbs model [1]_ estimates the diffuse fraction DF from global + horizontal irradiance through an empirical relationship between DF + and the ratio of GHI to extraterrestrial irradiance, Kt. The + function uses the diffuse fraction to compute DHI as + + .. math:: + + DHI = DF \times GHI + + DNI is then estimated as + + .. math:: + + DNI = (GHI - DHI)/\cos(Z) + + where Z is the zenith angle. + + Parameters + ---------- + ghi: numeric + Global horizontal irradiance in W/m^2. + zenith: numeric + True (not refraction-corrected) zenith angles in decimal degrees. + doy: scalar, array or DatetimeIndex + The day of the year. + + Returns + ------- + data : OrderedDict or DataFrame + Contains the following keys/columns: + + * ``dni``: the modeled direct normal irradiance in W/m^2. + * ``dhi``: the modeled diffuse horizontal irradiance in + W/m^2. + * ``kt``: Ratio of global to extraterrestrial irradiance + on a horizontal plane. + + References + ---------- + .. [1] D. G. Erbs, S. A. Klein and J. A. Duffie, Estimation of the + diffuse radiation fraction for hourly, daily and monthly-average + global radiation, Solar Energy 28(4), pp 293-302, 1982. Eq. 1 + + See also + -------- + dirint + disc + """ + + dni_extra = extraradiation(doy) + + # This Z needs to be the true Zenith angle, not apparent, + # to get extraterrestrial horizontal radiation) + i0_h = dni_extra * tools.cosd(zenith) + + kt = ghi / i0_h + kt = np.maximum(kt, 0) + + # For Kt <= 0.22, set the diffuse fraction + df = 1 - 0.09*kt + + # For Kt > 0.22 and Kt <= 0.8, set the diffuse fraction + df = np.where((kt > 0.22) & (kt <= 0.8), + 0.9511 - 0.1604*kt + 4.388*kt**2 - + 16.638*kt**3 + 12.336*kt**4, + df) + + # For Kt > 0.8, set the diffuse fraction + df = np.where(kt > 0.8, 0.165, df) + + dhi = df * ghi + + dni = (ghi - dhi) / tools.cosd(zenith) + + data = OrderedDict() + data['dni'] = dni + data['dhi'] = dhi + data['kt'] = kt + + if isinstance(dni, pd.Series): + data = pd.DataFrame(data) + + return data + + +def _get_perez_coefficients(perezmodel): + ''' + Find coefficients for the Perez model + + Parameters + ---------- + + perezmodel : string (optional, default='allsitescomposite1990') + + a character string which selects the desired set of Perez + coefficients. If model is not provided as an input, the default, + '1990' will be used. + + All possible model selections are: + + * '1990' + * 'allsitescomposite1990' (same as '1990') + * 'allsitescomposite1988' + * 'sandiacomposite1988' + * 'usacomposite1988' + * 'france1988' + * 'phoenix1988' + * 'elmonte1988' + * 'osage1988' + * 'albuquerque1988' + * 'capecanaveral1988' + * 'albany1988' + + Returns + -------- + F1coeffs, F2coeffs : (array, array) + F1 and F2 coefficients for the Perez model + + References + ---------- + [1] Loutzenhiser P.G. et. al. "Empirical validation of models to + compute solar irradiance on inclined surfaces for building energy + simulation" 2007, Solar Energy vol. 81. pp. 254-267 + + [2] Perez, R., Seals, R., Ineichen, P., Stewart, R., Menicucci, D., + 1987. A new simplified version of the Perez diffuse irradiance model + for tilted surfaces. Solar Energy 39(3), 221-232. + + [3] Perez, R., Ineichen, P., Seals, R., Michalsky, J., Stewart, R., + 1990. Modeling daylight availability and irradiance components from + direct and global irradiance. Solar Energy 44 (5), 271-289. + + [4] Perez, R. et. al 1988. "The Development and Verification of the + Perez Diffuse Radiation Model". SAND88-7030 + + ''' + coeffdict = { + 'allsitescomposite1990': [ + [-0.0080, 0.5880, -0.0620, -0.0600, 0.0720, -0.0220], + [0.1300, 0.6830, -0.1510, -0.0190, 0.0660, -0.0290], + [0.3300, 0.4870, -0.2210, 0.0550, -0.0640, -0.0260], + [0.5680, 0.1870, -0.2950, 0.1090, -0.1520, -0.0140], + [0.8730, -0.3920, -0.3620, 0.2260, -0.4620, 0.0010], + [1.1320, -1.2370, -0.4120, 0.2880, -0.8230, 0.0560], + [1.0600, -1.6000, -0.3590, 0.2640, -1.1270, 0.1310], + [0.6780, -0.3270, -0.2500, 0.1560, -1.3770, 0.2510]], + 'allsitescomposite1988': [ + [-0.0180, 0.7050, -0.071, -0.0580, 0.1020, -0.0260], + [0.1910, 0.6450, -0.1710, 0.0120, 0.0090, -0.0270], + [0.4400, 0.3780, -0.2560, 0.0870, -0.1040, -0.0250], + [0.7560, -0.1210, -0.3460, 0.1790, -0.3210, -0.0080], + [0.9960, -0.6450, -0.4050, 0.2600, -0.5900, 0.0170], + [1.0980, -1.2900, -0.3930, 0.2690, -0.8320, 0.0750], + [0.9730, -1.1350, -0.3780, 0.1240, -0.2580, 0.1490], + [0.6890, -0.4120, -0.2730, 0.1990, -1.6750, 0.2370]], + 'sandiacomposite1988': [ + [-0.1960, 1.0840, -0.0060, -0.1140, 0.1800, -0.0190], + [0.2360, 0.5190, -0.1800, -0.0110, 0.0200, -0.0380], + [0.4540, 0.3210, -0.2550, 0.0720, -0.0980, -0.0460], + [0.8660, -0.3810, -0.3750, 0.2030, -0.4030, -0.0490], + [1.0260, -0.7110, -0.4260, 0.2730, -0.6020, -0.0610], + [0.9780, -0.9860, -0.3500, 0.2800, -0.9150, -0.0240], + [0.7480, -0.9130, -0.2360, 0.1730, -1.0450, 0.0650], + [0.3180, -0.7570, 0.1030, 0.0620, -1.6980, 0.2360]], + 'usacomposite1988': [ + [-0.0340, 0.6710, -0.0590, -0.0590, 0.0860, -0.0280], + [0.2550, 0.4740, -0.1910, 0.0180, -0.0140, -0.0330], [0.4270, 0.3490, -0.2450, 0.0930, -0.1210, -0.0390], [0.7560, -0.2130, -0.3280, 0.1750, -0.3040, -0.0270], [1.0200, -0.8570, -0.3850, 0.2800, -0.6380, -0.0190], @@ -1294,274 +1576,12 @@ def _get_perez_coefficients(perezmodel): array = np.array(coeffdict[perezmodel]) - F1coeffs = array.T[0:3].T - F2coeffs = array.T[3:7].T + F1coeffs = array[:, 0:3] + F2coeffs = array[:, 3:7] return F1coeffs, F2coeffs -def disc(ghi, zenith, times, pressure=101325): - ''' - Estimate Direct Normal Irradiance from Global Horizontal Irradiance - using the DISC model. - - The DISC algorithm converts global horizontal irradiance to direct - normal irradiance through empirical relationships between the global - and direct clearness indices. - - Parameters - ---------- - - ghi : Series - Global horizontal irradiance in W/m^2. - - solar_zenith : Series - True (not refraction - corrected) solar zenith - angles in decimal degrees. - - times : DatetimeIndex - - pressure : float or Series - Site pressure in Pascal. - - Returns - ------- - DataFrame with the following keys: - * ``dni``: The modeled direct normal irradiance - in W/m^2 provided by the - Direct Insolation Simulation Code (DISC) model. - * ``kt``: Ratio of global to extraterrestrial - irradiance on a horizontal plane. - * ``airmass``: Airmass - - References - ---------- - - [1] Maxwell, E. L., "A Quasi-Physical Model for Converting Hourly - Global Horizontal to Direct Normal Insolation", Technical - Report No. SERI/TR-215-3087, Golden, CO: Solar Energy Research - Institute, 1987. - - [2] J.W. "Fourier series representation of the position of the sun". - Found at: - http://www.mail-archive.com/sundial@uni-koeln.de/msg01050.html on - January 12, 2012 - - See Also - -------- - atmosphere.alt2pres - dirint - ''' - - pvl_logger.debug('clearsky.disc') - - temp = pd.DataFrame(index=times, columns=['A','B','C'], dtype=float) - - doy = times.dayofyear - - DayAngle = 2. * np.pi*(doy - 1) / 365 - - re = (1.00011 + 0.034221*np.cos(DayAngle) + 0.00128*np.sin(DayAngle) - + 0.000719*np.cos(2.*DayAngle) + 7.7e-05*np.sin(2.*DayAngle) ) - - I0 = re * 1370. - I0h = I0 * np.cos(np.radians(zenith)) - - Ztemp = zenith.copy() - Ztemp[zenith > 87] = np.NaN - - AM = 1.0 / ( np.cos(np.radians(Ztemp)) + 0.15*( (93.885 - Ztemp)**(-1.253) ) ) * (pressure / 101325) - - Kt = ghi / I0h - Kt[Kt < 0] = 0 - Kt[Kt > 2] = np.NaN - - temp.A[Kt > 0.6] = -5.743 + 21.77*(Kt[Kt > 0.6]) - 27.49*(Kt[Kt > 0.6] ** 2) + 11.56*(Kt[Kt > 0.6] ** 3) - temp.B[Kt > 0.6] = 41.4 - 118.5*(Kt[Kt > 0.6]) + 66.05*(Kt[Kt > 0.6] ** 2) + 31.9*(Kt[Kt > 0.6] ** 3) - temp.C[Kt > 0.6] = -47.01 + 184.2*(Kt[Kt > 0.6]) - 222.0 * Kt[Kt > 0.6] ** 2 + 73.81*(Kt[Kt > 0.6] ** 3) - temp.A[Kt <= 0.6] = 0.512 - 1.56*(Kt[Kt <= 0.6]) + 2.286*(Kt[Kt <= 0.6] ** 2) - 2.222*(Kt[Kt <= 0.6] ** 3) - temp.B[Kt <= 0.6] = 0.37 + 0.962*(Kt[Kt <= 0.6]) - temp.C[Kt <= 0.6] = -0.28 + 0.932*(Kt[Kt <= 0.6]) - 2.048*(Kt[Kt <= 0.6] ** 2) - - delKn = temp.A + temp.B * np.exp(temp.C*AM) - - Knc = 0.866 - 0.122*(AM) + 0.0121*(AM ** 2) - 0.000653*(AM ** 3) + 1.4e-05*(AM ** 4) - Kn = Knc - delKn - - dni = Kn * I0 - - dni[zenith > 87] = np.NaN - dni[(ghi < 0) | (dni < 0)] = 0 - - dfout = pd.DataFrame({'dni':dni}) - dfout['kt'] = Kt - dfout['airmass'] = AM - - return dfout - - -def dirint(ghi, zenith, times, pressure=101325, use_delta_kt_prime=True, - temp_dew=None): - """ - Determine DNI from GHI using the DIRINT modification - of the DISC model. - - Implements the modified DISC model known as "DIRINT" introduced in [1]. - DIRINT predicts direct normal irradiance (DNI) from measured global - horizontal irradiance (GHI). DIRINT improves upon the DISC model by - using time-series GHI data and dew point temperature information. The - effectiveness of the DIRINT model improves with each piece of - information provided. - - Parameters - ---------- - ghi : pd.Series - Global horizontal irradiance in W/m^2. - - zenith : pd.Series - True (not refraction-corrected) zenith - angles in decimal degrees. If Z is a vector it must be of the - same size as all other vector inputs. Z must be >=0 and <=180. - - times : DatetimeIndex - - pressure : float or pd.Series - The site pressure in Pascal. - Pressure may be measured or an average pressure may be - calculated from site altitude. - - use_delta_kt_prime : bool - Indicates if the user would like to - utilize the time-series nature of the GHI measurements. A value of ``False`` - will not use the time-series improvements, any other numeric value - will use time-series improvements. It is recommended that time-series - data only be used if the time between measured data points is less - than 1.5 hours. If none of the input arguments are - vectors, then time-series improvements are not used (because it's not - a time-series). - - temp_dew : None, float, or pd.Series - Surface dew point temperatures, in degrees C. - Values of temp_dew may be numeric or NaN. Any - single time period point with a DewPtTemp=NaN does not have dew point - improvements applied. If DewPtTemp is not provided, then dew point - improvements are not applied. - - Returns - ------- - dni : pd.Series. - The modeled direct normal irradiance in W/m^2 provided by the - DIRINT model. - - References - ---------- - [1] Perez, R., P. Ineichen, E. Maxwell, R. Seals and A. Zelenka, (1992). - "Dynamic Global-to-Direct Irradiance Conversion Models". ASHRAE - Transactions-Research Series, pp. 354-369 - - [2] Maxwell, E. L., "A Quasi-Physical Model for Converting Hourly - Global Horizontal to Direct Normal Insolation", Technical - Report No. SERI/TR-215-3087, Golden, CO: Solar Energy Research - Institute, 1987. - - DIRINT model requires time series data (ie. one of the inputs must be a - vector of length >2. - """ - - pvl_logger.debug('clearsky.dirint') - - disc_out = disc(ghi, zenith, times, pressure=pressure) - kt = disc_out['kt'] - - # Absolute Airmass, per the DISC model - # Note that we calculate the AM pressure correction slightly differently - # than Perez. He uses altitude, we use pressure (which we calculate - # slightly differently) - airmass = (1./(tools.cosd(zenith) + 0.15*((93.885-zenith)**(-1.253))) * - pressure/101325) - - coeffs = _get_dirint_coeffs() - - kt_prime = kt / (1.031 * np.exp(-1.4/(0.9+9.4/airmass)) + 0.1) - kt_prime[kt_prime > 0.82] = 0.82 # From SRRL code. consider np.NaN - pvl_logger.debug('kt_prime:\n%s', kt_prime) - - # wholmgren: - # the use_delta_kt_prime statement is a port of the MATLAB code. - # I am confused by the abs() in the delta_kt_prime calculation. - # It is not the absolute value of the central difference. - if use_delta_kt_prime: - delta_kt_prime = 0.5*( (kt_prime - kt_prime.shift(1)).abs() - .add( - (kt_prime - kt_prime.shift(-1)).abs(), - fill_value=0)) - else: - delta_kt_prime = pd.Series(-1, index=times) - - if temp_dew is not None: - w = pd.Series(np.exp(0.07 * temp_dew - 0.075), index=times) - else: - w = pd.Series(-1, index=times) - - # @wholmgren: the following bin assignments use MATLAB's 1-indexing. - # Later, we'll subtract 1 to conform to Python's 0-indexing. - - # Create kt_prime bins - kt_prime_bin = pd.Series(0, index=times, dtype=np.int64) - kt_prime_bin[(kt_prime>=0) & (kt_prime<0.24)] = 1 - kt_prime_bin[(kt_prime>=0.24) & (kt_prime<0.4)] = 2 - kt_prime_bin[(kt_prime>=0.4) & (kt_prime<0.56)] = 3 - kt_prime_bin[(kt_prime>=0.56) & (kt_prime<0.7)] = 4 - kt_prime_bin[(kt_prime>=0.7) & (kt_prime<0.8)] = 5 - kt_prime_bin[(kt_prime>=0.8) & (kt_prime<=1)] = 6 - pvl_logger.debug('kt_prime_bin:\n%s', kt_prime_bin) - - # Create zenith angle bins - zenith_bin = pd.Series(0, index=times, dtype=np.int64) - zenith_bin[(zenith>=0) & (zenith<25)] = 1 - zenith_bin[(zenith>=25) & (zenith<40)] = 2 - zenith_bin[(zenith>=40) & (zenith<55)] = 3 - zenith_bin[(zenith>=55) & (zenith<70)] = 4 - zenith_bin[(zenith>=70) & (zenith<80)] = 5 - zenith_bin[(zenith>=80)] = 6 - pvl_logger.debug('zenith_bin:\n%s', zenith_bin) - - # Create the bins for w based on dew point temperature - w_bin = pd.Series(0, index=times, dtype=np.int64) - w_bin[(w>=0) & (w<1)] = 1 - w_bin[(w>=1) & (w<2)] = 2 - w_bin[(w>=2) & (w<3)] = 3 - w_bin[(w>=3)] = 4 - w_bin[(w == -1)] = 5 - pvl_logger.debug('w_bin:\n%s', w_bin) - - # Create delta_kt_prime binning. - delta_kt_prime_bin = pd.Series(0, index=times, dtype=np.int64) - delta_kt_prime_bin[(delta_kt_prime>=0) & (delta_kt_prime<0.015)] = 1 - delta_kt_prime_bin[(delta_kt_prime>=0.015) & (delta_kt_prime<0.035)] = 2 - delta_kt_prime_bin[(delta_kt_prime>=0.035) & (delta_kt_prime<0.07)] = 3 - delta_kt_prime_bin[(delta_kt_prime>=0.07) & (delta_kt_prime<0.15)] = 4 - delta_kt_prime_bin[(delta_kt_prime>=0.15) & (delta_kt_prime<0.3)] = 5 - delta_kt_prime_bin[(delta_kt_prime>=0.3) & (delta_kt_prime<=1)] = 6 - delta_kt_prime_bin[delta_kt_prime == -1] = 7 - pvl_logger.debug('delta_kt_prime_bin:\n%s', delta_kt_prime_bin) - - # subtract 1 to account for difference between MATLAB-style bin - # assignment and Python-style array lookup. - dirint_coeffs = coeffs[kt_prime_bin-1, zenith_bin-1, - delta_kt_prime_bin-1, w_bin-1] - # convert unassigned bins to nan - # use where to avoid boolean indexing deprecation - dirint_coeffs[np.where((kt_prime_bin == 0) | (zenith_bin == 0) | - (w_bin == 0) | (delta_kt_prime_bin == 0))] = np.nan - - dni = disc_out['dni'] * dirint_coeffs - - dni.name = 'dni' - - return dni - - def _get_dirint_coeffs(): """ A place to stash the dirint coefficients. @@ -1576,419 +1596,330 @@ def _get_dirint_coeffs(): # we create and assign values to an oversized array. # Then, we return the [1:, 1:, :, :] slice. - coeffs = np.zeros((7,7,7,5)) + coeffs = np.zeros((7, 7, 7, 5)) - coeffs[1,1,:,:] = [ + coeffs[1, 1, :, :] = [ [0.385230, 0.385230, 0.385230, 0.462880, 0.317440], [0.338390, 0.338390, 0.221270, 0.316730, 0.503650], [0.235680, 0.235680, 0.241280, 0.157830, 0.269440], [0.830130, 0.830130, 0.171970, 0.841070, 0.457370], [0.548010, 0.548010, 0.478000, 0.966880, 1.036370], [0.548010, 0.548010, 1.000000, 3.012370, 1.976540], - [0.582690, 0.582690, 0.229720, 0.892710, 0.569950 ]] + [0.582690, 0.582690, 0.229720, 0.892710, 0.569950]] - coeffs[1,2,:,:] = [ + coeffs[1, 2, :, :] = [ [0.131280, 0.131280, 0.385460, 0.511070, 0.127940], [0.223710, 0.223710, 0.193560, 0.304560, 0.193940], [0.229970, 0.229970, 0.275020, 0.312730, 0.244610], [0.090100, 0.184580, 0.260500, 0.687480, 0.579440], [0.131530, 0.131530, 0.370190, 1.380350, 1.052270], [1.116250, 1.116250, 0.928030, 3.525490, 2.316920], - [0.090100, 0.237000, 0.300040, 0.812470, 0.664970 ]] + [0.090100, 0.237000, 0.300040, 0.812470, 0.664970]] - coeffs[1,3,:,:] = [ + coeffs[1, 3, :, :] = [ [0.587510, 0.130000, 0.400000, 0.537210, 0.832490], [0.306210, 0.129830, 0.204460, 0.500000, 0.681640], [0.224020, 0.260620, 0.334080, 0.501040, 0.350470], [0.421540, 0.753970, 0.750660, 3.706840, 0.983790], [0.706680, 0.373530, 1.245670, 0.864860, 1.992630], [4.864400, 0.117390, 0.265180, 0.359180, 3.310820], - [0.392080, 0.493290, 0.651560, 1.932780, 0.898730 ]] + [0.392080, 0.493290, 0.651560, 1.932780, 0.898730]] - coeffs[1,4,:,:] = [ + coeffs[1, 4, :, :] = [ [0.126970, 0.126970, 0.126970, 0.126970, 0.126970], [0.810820, 0.810820, 0.810820, 0.810820, 0.810820], [3.241680, 2.500000, 2.291440, 2.291440, 2.291440], [4.000000, 3.000000, 2.000000, 0.975430, 1.965570], [12.494170, 12.494170, 8.000000, 5.083520, 8.792390], [21.744240, 21.744240, 21.744240, 21.744240, 21.744240], - [3.241680, 12.494170, 1.620760, 1.375250, 2.331620 ]] + [3.241680, 12.494170, 1.620760, 1.375250, 2.331620]] - coeffs[1,5,:,:] = [ + coeffs[1, 5, :, :] = [ [0.126970, 0.126970, 0.126970, 0.126970, 0.126970], [0.810820, 0.810820, 0.810820, 0.810820, 0.810820], [3.241680, 2.500000, 2.291440, 2.291440, 2.291440], [4.000000, 3.000000, 2.000000, 0.975430, 1.965570], [12.494170, 12.494170, 8.000000, 5.083520, 8.792390], [21.744240, 21.744240, 21.744240, 21.744240, 21.744240], - [3.241680, 12.494170, 1.620760, 1.375250, 2.331620 ]] + [3.241680, 12.494170, 1.620760, 1.375250, 2.331620]] - coeffs[1,6,:,:] = [ + coeffs[1, 6, :, :] = [ [0.126970, 0.126970, 0.126970, 0.126970, 0.126970], [0.810820, 0.810820, 0.810820, 0.810820, 0.810820], [3.241680, 2.500000, 2.291440, 2.291440, 2.291440], [4.000000, 3.000000, 2.000000, 0.975430, 1.965570], [12.494170, 12.494170, 8.000000, 5.083520, 8.792390], [21.744240, 21.744240, 21.744240, 21.744240, 21.744240], - [3.241680, 12.494170, 1.620760, 1.375250, 2.331620 ]] + [3.241680, 12.494170, 1.620760, 1.375250, 2.331620]] - coeffs[2,1,:,:] = [ + coeffs[2, 1, :, :] = [ [0.337440, 0.337440, 0.969110, 1.097190, 1.116080], [0.337440, 0.337440, 0.969110, 1.116030, 0.623900], [0.337440, 0.337440, 1.530590, 1.024420, 0.908480], [0.584040, 0.584040, 0.847250, 0.914940, 1.289300], [0.337440, 0.337440, 0.310240, 1.435020, 1.852830], [0.337440, 0.337440, 1.015010, 1.097190, 2.117230], - [0.337440, 0.337440, 0.969110, 1.145730, 1.476400 ]] + [0.337440, 0.337440, 0.969110, 1.145730, 1.476400]] - coeffs[2,2,:,:] = [ + coeffs[2, 2, :, :] = [ [0.300000, 0.300000, 0.700000, 1.100000, 0.796940], [0.219870, 0.219870, 0.526530, 0.809610, 0.649300], [0.386650, 0.386650, 0.119320, 0.576120, 0.685460], [0.746730, 0.399830, 0.470970, 0.986530, 0.785370], [0.575420, 0.936700, 1.649200, 1.495840, 1.335590], [1.319670, 4.002570, 1.276390, 2.644550, 2.518670], - [0.665190, 0.678910, 1.012360, 1.199940, 0.986580 ]] + [0.665190, 0.678910, 1.012360, 1.199940, 0.986580]] - coeffs[2,3,:,:] = [ + coeffs[2, 3, :, :] = [ [0.378870, 0.974060, 0.500000, 0.491880, 0.665290], [0.105210, 0.263470, 0.407040, 0.553460, 0.582590], [0.312900, 0.345240, 1.144180, 0.854790, 0.612280], [0.119070, 0.365120, 0.560520, 0.793720, 0.802600], [0.781610, 0.837390, 1.270420, 1.537980, 1.292950], [1.152290, 1.152290, 1.492080, 1.245370, 2.177100], - [0.424660, 0.529550, 0.966910, 1.033460, 0.958730 ]] + [0.424660, 0.529550, 0.966910, 1.033460, 0.958730]] - coeffs[2,4,:,:] = [ + coeffs[2, 4, :, :] = [ [0.310590, 0.714410, 0.252450, 0.500000, 0.607600], [0.975190, 0.363420, 0.500000, 0.400000, 0.502800], [0.175580, 0.196250, 0.476360, 1.072470, 0.490510], [0.719280, 0.698620, 0.657770, 1.190840, 0.681110], [0.426240, 1.464840, 0.678550, 1.157730, 0.978430], [2.501120, 1.789130, 1.387090, 2.394180, 2.394180], - [0.491640, 0.677610, 0.685610, 1.082400, 0.735410 ]] + [0.491640, 0.677610, 0.685610, 1.082400, 0.735410]] - coeffs[2,5,:,:] = [ + coeffs[2, 5, :, :] = [ [0.597000, 0.500000, 0.300000, 0.310050, 0.413510], [0.314790, 0.336310, 0.400000, 0.400000, 0.442460], [0.166510, 0.460440, 0.552570, 1.000000, 0.461610], [0.401020, 0.559110, 0.403630, 1.016710, 0.671490], [0.400360, 0.750830, 0.842640, 1.802600, 1.023830], [3.315300, 1.510380, 2.443650, 1.638820, 2.133990], - [0.530790, 0.745850, 0.693050, 1.458040, 0.804500 ]] + [0.530790, 0.745850, 0.693050, 1.458040, 0.804500]] - coeffs[2,6,:,:] = [ + coeffs[2, 6, :, :] = [ [0.597000, 0.500000, 0.300000, 0.310050, 0.800920], [0.314790, 0.336310, 0.400000, 0.400000, 0.237040], [0.166510, 0.460440, 0.552570, 1.000000, 0.581990], [0.401020, 0.559110, 0.403630, 1.016710, 0.898570], [0.400360, 0.750830, 0.842640, 1.802600, 3.400390], [3.315300, 1.510380, 2.443650, 1.638820, 2.508780], - [0.204340, 1.157740, 2.003080, 2.622080, 1.409380 ]] + [0.204340, 1.157740, 2.003080, 2.622080, 1.409380]] - coeffs[3,1,:,:] = [ + coeffs[3, 1, :, :] = [ [1.242210, 1.242210, 1.242210, 1.242210, 1.242210], [0.056980, 0.056980, 0.656990, 0.656990, 0.925160], [0.089090, 0.089090, 1.040430, 1.232480, 1.205300], [1.053850, 1.053850, 1.399690, 1.084640, 1.233340], [1.151540, 1.151540, 1.118290, 1.531640, 1.411840], [1.494980, 1.494980, 1.700000, 1.800810, 1.671600], - [1.018450, 1.018450, 1.153600, 1.321890, 1.294670 ]] + [1.018450, 1.018450, 1.153600, 1.321890, 1.294670]] - coeffs[3,2,:,:] = [ + coeffs[3, 2, :, :] = [ [0.700000, 0.700000, 1.023460, 0.700000, 0.945830], [0.886300, 0.886300, 1.333620, 0.800000, 1.066620], [0.902180, 0.902180, 0.954330, 1.126690, 1.097310], [1.095300, 1.075060, 1.176490, 1.139470, 1.096110], [1.201660, 1.201660, 1.438200, 1.256280, 1.198060], [1.525850, 1.525850, 1.869160, 1.985410, 1.911590], - [1.288220, 1.082810, 1.286370, 1.166170, 1.119330 ]] + [1.288220, 1.082810, 1.286370, 1.166170, 1.119330]] - coeffs[3,3,:,:] = [ + coeffs[3, 3, :, :] = [ [0.600000, 1.029910, 0.859890, 0.550000, 0.813600], [0.604450, 1.029910, 0.859890, 0.656700, 0.928840], [0.455850, 0.750580, 0.804930, 0.823000, 0.911000], [0.526580, 0.932310, 0.908620, 0.983520, 0.988090], [1.036110, 1.100690, 0.848380, 1.035270, 1.042380], [1.048440, 1.652720, 0.900000, 2.350410, 1.082950], - [0.817410, 0.976160, 0.861300, 0.974780, 1.004580 ]] + [0.817410, 0.976160, 0.861300, 0.974780, 1.004580]] - coeffs[3,4,:,:] = [ + coeffs[3, 4, :, :] = [ [0.782110, 0.564280, 0.600000, 0.600000, 0.665740], [0.894480, 0.680730, 0.541990, 0.800000, 0.669140], [0.487460, 0.818950, 0.841830, 0.872540, 0.709040], [0.709310, 0.872780, 0.908480, 0.953290, 0.844350], [0.863920, 0.947770, 0.876220, 1.078750, 0.936910], [1.280350, 0.866720, 0.769790, 1.078750, 0.975130], - [0.725420, 0.869970, 0.868810, 0.951190, 0.829220 ]] + [0.725420, 0.869970, 0.868810, 0.951190, 0.829220]] - coeffs[3,5,:,:] = [ + coeffs[3, 5, :, :] = [ [0.791750, 0.654040, 0.483170, 0.409000, 0.597180], [0.566140, 0.948990, 0.971820, 0.653570, 0.718550], [0.648710, 0.637730, 0.870510, 0.860600, 0.694300], [0.637630, 0.767610, 0.925670, 0.990310, 0.847670], [0.736380, 0.946060, 1.117590, 1.029340, 0.947020], [1.180970, 0.850000, 1.050000, 0.950000, 0.888580], - [0.700560, 0.801440, 0.961970, 0.906140, 0.823880 ]] + [0.700560, 0.801440, 0.961970, 0.906140, 0.823880]] - coeffs[3,6,:,:] = [ + coeffs[3, 6, :, :] = [ [0.500000, 0.500000, 0.586770, 0.470550, 0.629790], [0.500000, 0.500000, 1.056220, 1.260140, 0.658140], [0.500000, 0.500000, 0.631830, 0.842620, 0.582780], [0.554710, 0.734730, 0.985820, 0.915640, 0.898260], [0.712510, 1.205990, 0.909510, 1.078260, 0.885610], [1.899260, 1.559710, 1.000000, 1.150000, 1.120390], - [0.653880, 0.793120, 0.903320, 0.944070, 0.796130 ]] + [0.653880, 0.793120, 0.903320, 0.944070, 0.796130]] - coeffs[4,1,:,:] = [ + coeffs[4, 1, :, :] = [ [1.000000, 1.000000, 1.050000, 1.170380, 1.178090], [0.960580, 0.960580, 1.059530, 1.179030, 1.131690], [0.871470, 0.871470, 0.995860, 1.141910, 1.114600], [1.201590, 1.201590, 0.993610, 1.109380, 1.126320], [1.065010, 1.065010, 0.828660, 0.939970, 1.017930], [1.065010, 1.065010, 0.623690, 1.119620, 1.132260], - [1.071570, 1.071570, 0.958070, 1.114130, 1.127110 ]] + [1.071570, 1.071570, 0.958070, 1.114130, 1.127110]] - coeffs[4,2,:,:] = [ + coeffs[4, 2, :, :] = [ [0.950000, 0.973390, 0.852520, 1.092200, 1.096590], [0.804120, 0.913870, 0.980990, 1.094580, 1.042420], [0.737540, 0.935970, 0.999940, 1.056490, 1.050060], [1.032980, 1.034540, 0.968460, 1.032080, 1.015780], [0.900000, 0.977210, 0.945960, 1.008840, 0.969960], [0.600000, 0.750000, 0.750000, 0.844710, 0.899100], - [0.926800, 0.965030, 0.968520, 1.044910, 1.032310 ]] + [0.926800, 0.965030, 0.968520, 1.044910, 1.032310]] - coeffs[4,3,:,:] = [ + coeffs[4, 3, :, :] = [ [0.850000, 1.029710, 0.961100, 1.055670, 1.009700], [0.818530, 0.960010, 0.996450, 1.081970, 1.036470], [0.765380, 0.953500, 0.948260, 1.052110, 1.000140], [0.775610, 0.909610, 0.927800, 0.987800, 0.952100], [1.000990, 0.881880, 0.875950, 0.949100, 0.893690], [0.902370, 0.875960, 0.807990, 0.942410, 0.917920], - [0.856580, 0.928270, 0.946820, 1.032260, 0.972990 ]] + [0.856580, 0.928270, 0.946820, 1.032260, 0.972990]] - coeffs[4,4,:,:] = [ + coeffs[4, 4, :, :] = [ [0.750000, 0.857930, 0.983800, 1.056540, 0.980240], [0.750000, 0.987010, 1.013730, 1.133780, 1.038250], [0.800000, 0.947380, 1.012380, 1.091270, 0.999840], [0.800000, 0.914550, 0.908570, 0.999190, 0.915230], [0.778540, 0.800590, 0.799070, 0.902180, 0.851560], [0.680190, 0.317410, 0.507680, 0.388910, 0.646710], - [0.794920, 0.912780, 0.960830, 1.057110, 0.947950 ]] + [0.794920, 0.912780, 0.960830, 1.057110, 0.947950]] - coeffs[4,5,:,:] = [ + coeffs[4, 5, :, :] = [ [0.750000, 0.833890, 0.867530, 1.059890, 0.932840], [0.979700, 0.971470, 0.995510, 1.068490, 1.030150], [0.858850, 0.987920, 1.043220, 1.108700, 1.044900], [0.802400, 0.955110, 0.911660, 1.045070, 0.944470], [0.884890, 0.766210, 0.885390, 0.859070, 0.818190], [0.615680, 0.700000, 0.850000, 0.624620, 0.669300], - [0.835570, 0.946150, 0.977090, 1.049350, 0.979970 ]] + [0.835570, 0.946150, 0.977090, 1.049350, 0.979970]] - coeffs[4,6,:,:] = [ + coeffs[4, 6, :, :] = [ [0.689220, 0.809600, 0.900000, 0.789500, 0.853990], [0.854660, 0.852840, 0.938200, 0.923110, 0.955010], [0.938600, 0.932980, 1.010390, 1.043950, 1.041640], [0.843620, 0.981300, 0.951590, 0.946100, 0.966330], [0.694740, 0.814690, 0.572650, 0.400000, 0.726830], [0.211370, 0.671780, 0.416340, 0.297290, 0.498050], - [0.843540, 0.882330, 0.911760, 0.898420, 0.960210 ]] + [0.843540, 0.882330, 0.911760, 0.898420, 0.960210]] - coeffs[5,1,:,:] = [ + coeffs[5, 1, :, :] = [ [1.054880, 1.075210, 1.068460, 1.153370, 1.069220], [1.000000, 1.062220, 1.013470, 1.088170, 1.046200], [0.885090, 0.993530, 0.942590, 1.054990, 1.012740], [0.920000, 0.950000, 0.978720, 1.020280, 0.984440], [0.850000, 0.908500, 0.839940, 0.985570, 0.962180], [0.800000, 0.800000, 0.810080, 0.950000, 0.961550], - [1.038590, 1.063200, 1.034440, 1.112780, 1.037800 ]] + [1.038590, 1.063200, 1.034440, 1.112780, 1.037800]] - coeffs[5,2,:,:] = [ + coeffs[5, 2, :, :] = [ [1.017610, 1.028360, 1.058960, 1.133180, 1.045620], [0.920000, 0.998970, 1.033590, 1.089030, 1.022060], [0.912370, 0.949930, 0.979770, 1.020420, 0.981770], [0.847160, 0.935300, 0.930540, 0.955050, 0.946560], [0.880260, 0.867110, 0.874130, 0.972650, 0.883420], [0.627150, 0.627150, 0.700000, 0.774070, 0.845130], - [0.973700, 1.006240, 1.026190, 1.071960, 1.017240 ]] + [0.973700, 1.006240, 1.026190, 1.071960, 1.017240]] - coeffs[5,3,:,:] = [ + coeffs[5, 3, :, :] = [ [1.028710, 1.017570, 1.025900, 1.081790, 1.024240], [0.924980, 0.985500, 1.014100, 1.092210, 0.999610], [0.828570, 0.934920, 0.994950, 1.024590, 0.949710], [0.900810, 0.901330, 0.928830, 0.979570, 0.913100], [0.761030, 0.845150, 0.805360, 0.936790, 0.853460], [0.626400, 0.546750, 0.730500, 0.850000, 0.689050], - [0.957630, 0.985480, 0.991790, 1.050220, 0.987900 ]] + [0.957630, 0.985480, 0.991790, 1.050220, 0.987900]] - coeffs[5,4,:,:] = [ + coeffs[5, 4, :, :] = [ [0.992730, 0.993880, 1.017150, 1.059120, 1.017450], [0.975610, 0.987160, 1.026820, 1.075440, 1.007250], [0.871090, 0.933190, 0.974690, 0.979840, 0.952730], [0.828750, 0.868090, 0.834920, 0.905510, 0.871530], [0.781540, 0.782470, 0.767910, 0.764140, 0.795890], [0.743460, 0.693390, 0.514870, 0.630150, 0.715660], - [0.934760, 0.957870, 0.959640, 0.972510, 0.981640 ]] + [0.934760, 0.957870, 0.959640, 0.972510, 0.981640]] - coeffs[5,5,:,:] = [ + coeffs[5, 5, :, :] = [ [0.965840, 0.941240, 0.987100, 1.022540, 1.011160], [0.988630, 0.994770, 0.976590, 0.950000, 1.034840], [0.958200, 1.018080, 0.974480, 0.920000, 0.989870], [0.811720, 0.869090, 0.812020, 0.850000, 0.821050], [0.682030, 0.679480, 0.632450, 0.746580, 0.738550], [0.668290, 0.445860, 0.500000, 0.678920, 0.696510], - [0.926940, 0.953350, 0.959050, 0.876210, 0.991490 ]] + [0.926940, 0.953350, 0.959050, 0.876210, 0.991490]] - coeffs[5,6,:,:] = [ + coeffs[5, 6, :, :] = [ [0.948940, 0.997760, 0.850000, 0.826520, 0.998470], [1.017860, 0.970000, 0.850000, 0.700000, 0.988560], [1.000000, 0.950000, 0.850000, 0.606240, 0.947260], [1.000000, 0.746140, 0.751740, 0.598390, 0.725230], [0.922210, 0.500000, 0.376800, 0.517110, 0.548630], [0.500000, 0.450000, 0.429970, 0.404490, 0.539940], - [0.960430, 0.881630, 0.775640, 0.596350, 0.937680 ]] + [0.960430, 0.881630, 0.775640, 0.596350, 0.937680]] - coeffs[6,1,:,:] = [ + coeffs[6, 1, :, :] = [ [1.030000, 1.040000, 1.000000, 1.000000, 1.049510], [1.050000, 0.990000, 0.990000, 0.950000, 0.996530], [1.050000, 0.990000, 0.990000, 0.820000, 0.971940], [1.050000, 0.790000, 0.880000, 0.820000, 0.951840], [1.000000, 0.530000, 0.440000, 0.710000, 0.928730], [0.540000, 0.470000, 0.500000, 0.550000, 0.773950], - [1.038270, 0.920180, 0.910930, 0.821140, 1.034560 ]] + [1.038270, 0.920180, 0.910930, 0.821140, 1.034560]] - coeffs[6,2,:,:] = [ + coeffs[6, 2, :, :] = [ [1.041020, 0.997520, 0.961600, 1.000000, 1.035780], [0.948030, 0.980000, 0.900000, 0.950360, 0.977460], [0.950000, 0.977250, 0.869270, 0.800000, 0.951680], [0.951870, 0.850000, 0.748770, 0.700000, 0.883850], [0.900000, 0.823190, 0.727450, 0.600000, 0.839870], [0.850000, 0.805020, 0.692310, 0.500000, 0.788410], - [1.010090, 0.895270, 0.773030, 0.816280, 1.011680 ]] + [1.010090, 0.895270, 0.773030, 0.816280, 1.011680]] - coeffs[6,3,:,:] = [ + coeffs[6, 3, :, :] = [ [1.022450, 1.004600, 0.983650, 1.000000, 1.032940], [0.943960, 0.999240, 0.983920, 0.905990, 0.978150], [0.936240, 0.946480, 0.850000, 0.850000, 0.930320], [0.816420, 0.885000, 0.644950, 0.817650, 0.865310], [0.742960, 0.765690, 0.561520, 0.700000, 0.827140], [0.643870, 0.596710, 0.474460, 0.600000, 0.651200], - [0.971740, 0.940560, 0.714880, 0.864380, 1.001650 ]] + [0.971740, 0.940560, 0.714880, 0.864380, 1.001650]] - coeffs[6,4,:,:] = [ + coeffs[6, 4, :, :] = [ [0.995260, 0.977010, 1.000000, 1.000000, 1.035250], [0.939810, 0.975250, 0.939980, 0.950000, 0.982550], [0.876870, 0.879440, 0.850000, 0.900000, 0.917810], [0.873480, 0.873450, 0.751470, 0.850000, 0.863040], [0.761470, 0.702360, 0.638770, 0.750000, 0.783120], [0.734080, 0.650000, 0.600000, 0.650000, 0.715660], - [0.942160, 0.919100, 0.770340, 0.731170, 0.995180 ]] + [0.942160, 0.919100, 0.770340, 0.731170, 0.995180]] - coeffs[6,5,:,:] = [ + coeffs[6, 5, :, :] = [ [0.952560, 0.916780, 0.920000, 0.900000, 1.005880], [0.928620, 0.994420, 0.900000, 0.900000, 0.983720], [0.913070, 0.850000, 0.850000, 0.800000, 0.924280], [0.868090, 0.807170, 0.823550, 0.600000, 0.844520], [0.769570, 0.719870, 0.650000, 0.550000, 0.733500], [0.580250, 0.650000, 0.600000, 0.500000, 0.628850], - [0.904770, 0.852650, 0.708370, 0.493730, 0.949030 ]] + [0.904770, 0.852650, 0.708370, 0.493730, 0.949030]] - coeffs[6,6,:,:] = [ + coeffs[6, 6, :, :] = [ [0.911970, 0.800000, 0.800000, 0.800000, 0.956320], [0.912620, 0.682610, 0.750000, 0.700000, 0.950110], [0.653450, 0.659330, 0.700000, 0.600000, 0.856110], [0.648440, 0.600000, 0.641120, 0.500000, 0.695780], [0.570000, 0.550000, 0.598800, 0.400000, 0.560150], [0.475230, 0.500000, 0.518640, 0.339970, 0.520230], - [0.743440, 0.592190, 0.603060, 0.316930, 0.794390 ]] - - return coeffs[1:,1:,:,:] - - -def erbs(ghi, zenith, doy): - r""" - Estimate DNI and DHI from GHI using the Erbs model. - - The Erbs model [1]_ estimates the diffuse fraction DF from global - horizontal irradiance through an empirical relationship between DF - and the ratio of GHI to extraterrestrial irradiance, Kt. The - function uses the diffuse fraction to compute DHI as - - .. math:: - - DHI = DF \times GHI - - DNI is then estimated as - - .. math:: - - DNI = (GHI - DHI)/\cos(Z) - - where Z is the zenith angle. + [0.743440, 0.592190, 0.603060, 0.316930, 0.794390]] - Parameters - ---------- - ghi: scalar, array or Series - Global horizontal irradiance in W/m^2. - zenith: scalar, array or Series - True (not refraction-corrected) zenith angles in decimal degrees. - doy: scalar, array or DatetimeIndex - The day of the year. - - Returns - ------- - data : DataFrame - The DataFrame will have the following columns: - - * ``dni``: the modeled direct normal irradiance in W/m^2. - * ``dhi``: the modeled diffuse horizontal irradiance in W/m^2. - * ``kt``: Ratio of global to extraterrestrial irradiance - on a horizontal plane. - - References - ---------- - .. [1] D. G. Erbs, S. A. Klein and J. A. Duffie, Estimation of the - diffuse radiation fraction for hourly, daily and monthly-average - global radiation, Solar Energy 28(4), pp 293-302, 1982. Eq. 1 - - See also - -------- - dirint - disc - """ - - # enables all scalar input - ghi = pd.Series(ghi) - zenith = pd.Series(zenith) - - dni_extra = extraradiation(doy) - - # This Z needs to be the true Zenith angle, not apparent, - # to get extraterrestrial horizontal radiation) - i0_h = dni_extra * tools.cosd(zenith) - - kt = ghi / i0_h - kt[kt < 0] = 0 - - # For Kt <= 0.22, set the diffuse fraction - df = 1 - 0.09*kt - - # For Kt > 0.22 and Kt <= 0.8, set the diffuse fraction - mask = (kt > 0.22) & (kt <= 0.8) - df[mask] = ( - 0.9511 - 0.1604*kt[mask] + 4.388*kt[mask]**2 - 16.638*kt[mask]**3 + - 12.336*kt[mask]**4) - - # For Kt > 0.8, set the diffuse fraction - df[kt > 0.8] = 0.165 - - dhi = df * ghi - - dni = (ghi - dhi) / tools.cosd(zenith) - - data = pd.DataFrame() - data['dni'] = dni - data['dhi'] = dhi - data['kt'] = kt - - return data + return coeffs[1:, 1:, :, :] diff --git a/pvlib/test/conftest.py b/pvlib/test/conftest.py index 412adfeb97..a128f17091 100644 --- a/pvlib/test/conftest.py +++ b/pvlib/test/conftest.py @@ -24,8 +24,15 @@ requires_ephem = pytest.mark.skipif(not has_ephem, reason='requires ephem') -incompatible_pandas_0131 = pytest.mark.skipif( - pd.__version__ == '0.13.1', reason='requires numpy 1.10 or greater') +def pandas_0_17(): + version = tuple(map(int, pd.__version__.split('.'))) + if version[0] <= 0 and version[1] < 17: + return False + else: + return True + +needs_pandas_0_17 = pytest.mark.skipif( + not pandas_0_17(), reason='requires pandas 0.17 or greater') def numpy_1_10(): diff --git a/pvlib/test/test_atmosphere.py b/pvlib/test/test_atmosphere.py index d06ecd507b..5522fa337b 100644 --- a/pvlib/test/test_atmosphere.py +++ b/pvlib/test/test_atmosphere.py @@ -12,16 +12,11 @@ from pvlib import atmosphere -# setup times and location to be tested. -times = pd.date_range(start=datetime.datetime(2014,6,24), - end=datetime.datetime(2014,6,26), freq='1Min') +latitude, longitude, tz, altitude = 32.2, -111, 'US/Arizona', 700 -tus = Location(32.2, -111, 'US/Arizona', 700) +times = pd.date_range(start='20140626', end='20140626', freq='6h', tz=tz) -times_localized = times.tz_localize(tus.tz) - -ephem_data = solarposition.get_solarposition(times_localized, tus.latitude, - tus.longitude) +ephem_data = solarposition.get_solarposition(times, latitude, longitude) # need to add physical tests instead of just functional tests @@ -38,7 +33,18 @@ def test_alt2press(): ['simple', 'kasten1966', 'youngirvine1967', 'kastenyoung1989', 'gueymard1993', 'young1994', 'pickering2002']) def test_airmass(model): - atmosphere.relativeairmass(ephem_data['zenith'], model) + out = atmosphere.relativeairmass(ephem_data['zenith'], model) + assert isinstance(out, pd.Series) + out = atmosphere.relativeairmass(ephem_data['zenith'].values, model) + assert isinstance(out, np.ndarray) + + +def test_airmass_scalar(): + assert not np.isnan(atmosphere.relativeairmass(10)) + + +def test_airmass_scalar_nan(): + assert np.isnan(atmosphere.relativeairmass(100)) def test_airmass_invalid(): diff --git a/pvlib/test/test_irradiance.py b/pvlib/test/test_irradiance.py index 67d3c24e5d..3b8b8e4f51 100644 --- a/pvlib/test/test_irradiance.py +++ b/pvlib/test/test_irradiance.py @@ -1,10 +1,12 @@ +from collections import OrderedDict + import numpy as np import pandas as pd import pytest from numpy.testing import assert_almost_equal, assert_allclose -from pandas.util.testing import assert_frame_equal +from pandas.util.testing import assert_frame_equal, assert_series_equal from pvlib.location import Location from pvlib import clearsky @@ -18,8 +20,7 @@ tus = Location(32.2, -111, 'US/Arizona', 700) # must include night values -times = pd.date_range(start='20140624', end='20140626', freq='1Min', - tz=tus.tz) +times = pd.date_range(start='20140624', freq='6H', periods=4, tz=tus.tz) ephem_data = solarposition.get_solarposition( times, tus.latitude, tus.longitude, method='nrel_numpy') @@ -135,10 +136,29 @@ def test_king(): def test_perez(): am = atmosphere.relativeairmass(ephem_data['apparent_zenith']) - out = irradiance.perez(40, 180, irrad_data['dhi'], irrad_data['dni'], + dni = irrad_data['dni'].copy() + dni.iloc[2] = np.nan + out = irradiance.perez(40, 180, irrad_data['dhi'], dni, dni_et, ephem_data['apparent_zenith'], ephem_data['azimuth'], am) - assert not out.isnull().any() + expected = pd.Series(np.array( + [ 0. , 31.46046871, np.nan, 45.45539877]), + index=times) + assert_series_equal(out, expected) + + + +def test_perez_arrays(): + am = atmosphere.relativeairmass(ephem_data['apparent_zenith']) + dni = irrad_data['dni'].copy() + dni.iloc[2] = np.nan + out = irradiance.perez(40, 180, irrad_data['dhi'].values, dni.values, + dni_et, ephem_data['apparent_zenith'].values, + ephem_data['azimuth'].values, am.values) + expected = np.array( + [ 0. , 31.46046871, np.nan, 45.45539877]) + assert_allclose(out, expected) + # klutcher (misspelling) will be removed in 0.3 def test_total_irrad(): @@ -161,6 +181,25 @@ def test_total_irrad(): 'poa_ground_diffuse'] +@pytest.mark.parametrize('model', ['isotropic', 'klucher', + 'haydavies', 'reindl', 'king', 'perez']) +def test_total_irrad_scalars(model): + total = irradiance.total_irrad( + 32, 180, + 10, 180, + dni=1000, ghi=1100, + dhi=100, + dni_extra=1400, airmass=1, + model=model, + surface_type='urban') + + assert list(total.keys()) == ['poa_global', 'poa_direct', + 'poa_diffuse', 'poa_sky_diffuse', + 'poa_ground_diffuse'] + # test that none of the values are nan + assert np.isnan(np.array(list(total.values()))).sum() == 0 + + def test_globalinplane(): aoi = irradiance.aoi(40, 180, ephem_data['apparent_zenith'], ephem_data['azimuth']) @@ -270,10 +309,13 @@ def test_erbs_all_scalar(): ghi = 1000 zenith = 10 doy = 180 - expected = pd.DataFrame(np. - array([[ 8.42358014e+02, 1.70439297e+02, 7.68919470e-01]]), - columns=['dni', 'dhi', 'kt']) + + expected = OrderedDict() + expected['dni'] = 8.42358014e+02 + expected['dhi'] = 1.70439297e+02 + expected['kt'] = 7.68919470e-01 out = irradiance.erbs(ghi, zenith, doy) - assert_frame_equal(out, expected) + for k, v in out.items(): + assert_allclose(v, expected[k], 5) diff --git a/pvlib/test/test_location.py b/pvlib/test/test_location.py index 7be7543861..e0bcea859e 100644 --- a/pvlib/test/test_location.py +++ b/pvlib/test/test_location.py @@ -12,6 +12,7 @@ from pvlib.location import Location from test_solarposition import expected_solpos +from conftest import requires_scipy aztz = pytz.timezone('US/Arizona') @@ -50,6 +51,7 @@ def test_location_print_pytz(): assert tus.__str__() == expected_str +@requires_scipy def test_get_clearsky(): tus = Location(32.2, -111, 'US/Arizona', 700, 'Tucson') times = pd.DatetimeIndex(start='20160101T0600-0700', diff --git a/pvlib/test/test_modelchain.py b/pvlib/test/test_modelchain.py index 1d166aade4..7168fa4bb3 100644 --- a/pvlib/test/test_modelchain.py +++ b/pvlib/test/test_modelchain.py @@ -12,6 +12,7 @@ import pytest from test_pvsystem import sam_data +from conftest import requires_scipy @pytest.fixture @@ -58,6 +59,7 @@ def run_orientation_strategy(strategy, expected): assert system.surface_azimuth == expected[1] +@requires_scipy def test_run_model(system, location): mc = ModelChain(system, location) times = pd.date_range('20160101 1200-0700', periods=2, freq='6H') @@ -105,6 +107,7 @@ def test_run_model_gueymard_perez(system, location): assert_series_equal(ac, expected) +@requires_scipy def test_run_model_with_weather(system, location): mc = ModelChain(system, location) times = pd.date_range('20160101 1200-0700', periods=2, freq='6H') @@ -116,6 +119,7 @@ def test_run_model_with_weather(system, location): assert_series_equal(ac, expected) +@requires_scipy def test_run_model_tracker(system, location): system = SingleAxisTracker(module_parameters=system.module_parameters, inverter_parameters=system.inverter_parameters) @@ -140,6 +144,7 @@ def test_bad_get_orientation(): modelchain.get_orientation('bad value') +@requires_scipy def test_basic_chain_required(sam_data): times = pd.DatetimeIndex(start='20160101 1200-0700', end='20160101 1800-0700', freq='6H') @@ -156,6 +161,7 @@ def test_basic_chain_required(sam_data): altitude=altitude) +@requires_scipy def test_basic_chain_alt_az(sam_data): times = pd.DatetimeIndex(start='20160101 1200-0700', end='20160101 1800-0700', freq='6H') @@ -179,6 +185,7 @@ def test_basic_chain_alt_az(sam_data): assert_series_equal(ac, expected) +@requires_scipy def test_basic_chain_strategy(sam_data): times = pd.DatetimeIndex(start='20160101 1200-0700', end='20160101 1800-0700', freq='6H') @@ -200,6 +207,7 @@ def test_basic_chain_strategy(sam_data): assert_series_equal(ac, expected) +@requires_scipy def test_basic_chain_altitude_pressure(sam_data): times = pd.DatetimeIndex(start='20160101 1200-0700', end='20160101 1800-0700', freq='6H') diff --git a/pvlib/test/test_pvsystem.py b/pvlib/test/test_pvsystem.py index 100e0c0274..dda74d0f17 100644 --- a/pvlib/test/test_pvsystem.py +++ b/pvlib/test/test_pvsystem.py @@ -18,7 +18,7 @@ from pvlib import solarposition from pvlib.location import Location -from conftest import needs_numpy_1_10 +from conftest import needs_numpy_1_10, requires_scipy latitude = 32.2 longitude = -111 @@ -234,21 +234,25 @@ def test_PVSystem_calcparams_desoto(sam_data): assert_allclose(nNsVth, 0.473) +@requires_scipy def test_v_from_i(): output = pvsystem.v_from_i(20, .1, .5, 3, 6e-7, 7) assert_allclose(7.5049875193450521, output, 5) +@requires_scipy def test_v_from_i_big(): output = pvsystem.v_from_i(500, 10, 4.06, 0, 6e-10, 1.2) assert_allclose(86.320000493521079, output, 5) +@requires_scipy def test_i_from_v(): output = pvsystem.i_from_v(20, .1, .5, 40, 6e-7, 7) assert_allclose(-299.746389916, output, 5) +@requires_scipy def test_PVSystem_i_from_v(sam_data): module = 'Example_Module' module_parameters = sam_data['cecmod'][module] @@ -258,6 +262,7 @@ def test_PVSystem_i_from_v(sam_data): assert_allclose(-299.746389916, output, 5) +@requires_scipy def test_singlediode_series(sam_data): module = 'Example_Module' module_parameters = sam_data['cecmod'][module] @@ -274,6 +279,7 @@ def test_singlediode_series(sam_data): assert isinstance(out, pd.DataFrame) +@requires_scipy def test_singlediode_floats(sam_data): module = 'Example_Module' module_parameters = sam_data['cecmod'][module] @@ -290,6 +296,7 @@ def test_singlediode_floats(sam_data): assert_allclose(expected[k], v, atol=3) +@requires_scipy def test_PVSystem_singlediode_floats(sam_data): module = 'Example_Module' module_parameters = sam_data['cecmod'][module] diff --git a/pvlib/test/test_solarposition.py b/pvlib/test/test_solarposition.py index 3bb46dc80d..fae0da1534 100644 --- a/pvlib/test/test_solarposition.py +++ b/pvlib/test/test_solarposition.py @@ -10,7 +10,7 @@ from pvlib.location import Location from pvlib import solarposition -from conftest import (requires_ephem, incompatible_pandas_0131, +from conftest import (requires_ephem, needs_pandas_0_17, requires_spa_c, requires_numba) @@ -116,7 +116,7 @@ def test_spa_python_numba_physical_dst(expected_solpos): assert_frame_equal(expected_solpos, ephem_data[expected_solpos.columns]) -@incompatible_pandas_0131 +@needs_pandas_0_17 def test_get_sun_rise_set_transit(): south = Location(-35.0, 0.0, tz='UTC') times = pd.DatetimeIndex([datetime.datetime(1996, 7, 5, 0), diff --git a/setup.py b/setup.py index 822160bb49..3709815d80 100755 --- a/setup.py +++ b/setup.py @@ -40,8 +40,8 @@ MAINTAINER_EMAIL = 'holmgren@email.arizona.edu' URL = 'https://github.com/pvlib/pvlib-python' -INSTALL_REQUIRES = ['numpy >= 1.8.2', - 'pandas >= 0.13.1', +INSTALL_REQUIRES = ['numpy >= 1.9.0', + 'pandas >= 0.14.0', 'pytz', 'six', ]