From f773e2dc748fc134a736570e1b8876de1021deb3 Mon Sep 17 00:00:00 2001 From: Will Holmgren Date: Wed, 15 Jun 2016 11:34:45 -0700 Subject: [PATCH] implement simplified solis clear sky model (#148) * implement simplified solis. needs tests * fix index stripping * add tests. coerce pw to min 0.2 * add simplified solis option to location.get_clearsky. reorder simplified solis output to match ineichen * update whatsnew * remove return_raw in favor of OrderedDict. improve nan handling --- docs/sphinx/source/whatsnew/v0.3.3.txt | 3 +- pvlib/clearsky.py | 290 +++++++++++++++++++++---- pvlib/location.py | 133 +++++++----- pvlib/solarposition.py | 2 +- pvlib/test/test_clearsky.py | 218 ++++++++++++++++++- pvlib/test/test_location.py | 101 ++++++++- 6 files changed, 640 insertions(+), 107 deletions(-) diff --git a/docs/sphinx/source/whatsnew/v0.3.3.txt b/docs/sphinx/source/whatsnew/v0.3.3.txt index a6bc3d61e2..47409e9dfb 100644 --- a/docs/sphinx/source/whatsnew/v0.3.3.txt +++ b/docs/sphinx/source/whatsnew/v0.3.3.txt @@ -1,6 +1,6 @@ .. _whatsnew_0330: -v0.3.3 (June xx, 2016) +v0.3.3 (June 15, 2016) ----------------------- This is a minor release from 0.3.2. @@ -32,6 +32,7 @@ Enhancements (:issue:`190`) * Improve speed of ``singlediode`` function by using ``v_from_i`` to determine ``v_oc``. Speed is ~2x faster. (:issue:`190`) +* Adds the Simplified Solis clear sky model. (:issue:`148`) Bug fixes diff --git a/pvlib/clearsky.py b/pvlib/clearsky.py index 7a0051679f..e227b6ee92 100644 --- a/pvlib/clearsky.py +++ b/pvlib/clearsky.py @@ -1,5 +1,5 @@ """ -The ``clearsky`` module contains several methods +The ``clearsky`` module contains several methods to calculate clear sky GHI, DNI, and DHI. """ @@ -9,6 +9,7 @@ logger = logging.getLogger('pvlib') import os +from collections import OrderedDict import numpy as np import pandas as pd @@ -19,7 +20,7 @@ from pvlib import solarposition -def ineichen(time, latitude, longitude, altitude=0, linke_turbidity=None, +def ineichen(time, latitude, longitude, altitude=0, linke_turbidity=None, solarposition_method='nrel_numpy', zenith_data=None, airmass_model='young1994', airmass_data=None, interp_turbidity=True): @@ -31,37 +32,37 @@ def ineichen(time, latitude, longitude, altitude=0, linke_turbidity=None, the clear-sky diffuse horizontal (DHI) component as the difference between GHI and DNI*cos(zenith) as presented in [1, 2]. A report on clear sky models found the Ineichen/Perez model to have excellent performance - with a minimal input data set [3]. - + with a minimal input data set [3]. + Default values for montly Linke turbidity provided by SoDa [4, 5]. Parameters ----------- time : pandas.DatetimeIndex - + latitude : float - + longitude : float altitude : float - + linke_turbidity : None or float If None, uses ``LinkeTurbidities.mat`` lookup table. - + solarposition_method : string - Sets the solar position algorithm. + Sets the solar position algorithm. See solarposition.get_solarposition() - + zenith_data : None or Series If None, ephemeris data will be calculated using ``solarposition_method``. - + airmass_model : string See pvlib.airmass.relativeairmass(). - + airmass_data : None or Series - If None, absolute air mass data will be calculated using + If None, absolute air mass data will be calculated using ``airmass_model`` and location.alitude. - + interp_turbidity : bool If ``True``, interpolates the monthly Linke turbidity values found in ``LinkeTurbidities.mat`` to daily values. @@ -100,16 +101,16 @@ def ineichen(time, latitude, longitude, altitude=0, linke_turbidity=None, # Initial implementation of this algorithm by Matthew Reno. # Ported to python by Rob Andrews # Added functionality by Will Holmgren (@wholmgren) - + I0 = irradiance.extraradiation(time.dayofyear) - + if zenith_data is None: ephem_data = solarposition.get_solarposition(time, latitude=latitude, longitude=longitude, altitude=altitude, method=solarposition_method) - time = ephem_data.index # fixes issue with time possibly not being tz-aware + time = ephem_data.index # fixes issue with time possibly not being tz-aware try: ApparentZenith = ephem_data['apparent_zenith'] except KeyError: @@ -119,7 +120,7 @@ def ineichen(time, latitude, longitude, altitude=0, linke_turbidity=None, ApparentZenith = zenith_data #ApparentZenith[ApparentZenith >= 90] = 90 # can cause problems in edge cases - + if linke_turbidity is None: TL = lookup_linke_turbidity(time, latitude, longitude, interp_turbidity=interp_turbidity) @@ -128,13 +129,13 @@ def ineichen(time, latitude, longitude, altitude=0, linke_turbidity=None, # Get the absolute airmass assuming standard local pressure (per # alt2pres) using Kasten and Young's 1989 formula for airmass. - + if airmass_data is None: AMabsolute = atmosphere.absoluteairmass(airmass_relative=atmosphere.relativeairmass(ApparentZenith, airmass_model), pressure=atmosphere.alt2pres(altitude)) else: AMabsolute = airmass_data - + fh1 = np.exp(-altitude/8000.) fh2 = np.exp(-altitude/1250.) cg1 = 5.09e-05 * altitude + 0.868 @@ -156,37 +157,37 @@ def ineichen(time, latitude, longitude, altitude=0, linke_turbidity=None, # TLcorr = TL; # TLcorr(TL < 2) = TLcorr(TL < 2) - 0.25 .* (2-TLcorr(TL < 2)) .^ (0.5); - # This equation is found in Solar Energy 73, pg 311. + # This equation is found in Solar Energy 73, pg 311. # Full ref: Perez et. al., Vol. 73, pp. 307-317 (2002). - # It is slightly different than the equation given in Solar Energy 73, pg 156. - # We used the equation from pg 311 because of the existence of known typos - # in the pg 156 publication (notably the fh2-(TL-1) should be fh2 * (TL-1)). - + # It is slightly different than the equation given in Solar Energy 73, pg 156. + # We used the equation from pg 311 because of the existence of known typos + # in the pg 156 publication (notably the fh2-(TL-1) should be fh2 * (TL-1)). + cos_zenith = tools.cosd(ApparentZenith) - + clearsky_GHI = ( cg1 * I0 * cos_zenith * np.exp(-cg2*AMabsolute*(fh1 + fh2*(TL - 1))) * np.exp(0.01*AMabsolute**1.8) ) clearsky_GHI[clearsky_GHI < 0] = 0 - + # BncI == "normal beam clear sky radiation" b = 0.664 + 0.163/fh1 BncI = b * I0 * np.exp( -0.09 * AMabsolute * (TL - 1) ) logger.debug('b=%s', b) - + # "empirical correction" SE 73, 157 & SE 73, 312. BncI_2 = ( clearsky_GHI * ( 1 - (0.1 - 0.2*np.exp(-TL))/(0.1 + 0.882/fh1) ) / cos_zenith ) clearsky_DNI = np.minimum(BncI, BncI_2) - + clearsky_DHI = clearsky_GHI - clearsky_DNI*cos_zenith - - df_out = pd.DataFrame({'ghi':clearsky_GHI, 'dni':clearsky_DNI, + + df_out = pd.DataFrame({'ghi':clearsky_GHI, 'dni':clearsky_DNI, 'dhi':clearsky_DHI}) df_out.fillna(0, inplace=True) - + return df_out @@ -222,11 +223,11 @@ def lookup_linke_turbidity(time, latitude, longitude, filepath=None, # from -180 to 180; and the depth (third dimension) represents months of # the year from January (1) to December (12). To determine the Linke # turbidity for a position on the Earth's surface for a given month do the - # following: LT = LinkeTurbidity(LatitudeIndex, LongitudeIndex, month). - # Note that the numbers within the matrix are 20 * Linke Turbidity, + # following: LT = LinkeTurbidity(LatitudeIndex, LongitudeIndex, month). + # Note that the numbers within the matrix are 20 * Linke Turbidity, # so divide the number from the file by 20 to get the - # turbidity. - + # turbidity. + try: import scipy.io except ImportError: @@ -252,7 +253,7 @@ def lookup_linke_turbidity(time, latitude, longitude, filepath=None, # Assume that data corresponds to the value at # the middle of each month. # This means that we need to add previous Dec and next Jan - # to the array so that the interpolation will work for + # to the array so that the interpolation will work for # Jan 1 - Jan 15 and Dec 16 - Dec 31. # Then we map the month value to the day of year value. # This is approximate and could be made more accurate. @@ -274,12 +275,12 @@ def lookup_linke_turbidity(time, latitude, longitude, filepath=None, def haurwitz(apparent_zenith): ''' Determine clear sky GHI from Haurwitz model. - + Implements the Haurwitz clear sky model for global horizontal irradiance (GHI) as presented in [1, 2]. A report on clear sky models found the Haurwitz model to have the best performance of models which require only zenith angle [3]. Extreme care should - be taken in the interpretation of this result! + be taken in the interpretation of this result! Parameters ---------- @@ -288,7 +289,7 @@ def haurwitz(apparent_zenith): in degrees. Returns - ------- + ------- pd.Series The modeled global horizonal irradiance in W/m^2 provided by the Haurwitz clear-sky model. @@ -298,10 +299,10 @@ def haurwitz(apparent_zenith): References ---------- - [1] B. Haurwitz, "Insolation in Relation to Cloudiness and Cloud + [1] B. Haurwitz, "Insolation in Relation to Cloudiness and Cloud Density," Journal of Meteorology, vol. 2, pp. 154-166, 1945. - [2] B. Haurwitz, "Insolation in Relation to Cloud Type," Journal of + [2] B. Haurwitz, "Insolation in Relation to Cloud Type," Journal of Meteorology, vol. 3, pp. 123-124, 1946. [3] M. Reno, C. Hansen, and J. Stein, "Global Horizontal Irradiance Clear @@ -314,16 +315,215 @@ def haurwitz(apparent_zenith): clearsky_GHI = 1098.0 * cos_zenith * np.exp(-0.059/cos_zenith) clearsky_GHI[clearsky_GHI < 0] = 0 - + df_out = pd.DataFrame({'ghi':clearsky_GHI}) - + return df_out def _linearly_scale(inputmatrix, inputmin, inputmax, outputmin, outputmax): """ used by linke turbidity lookup function """ - + inputrange = inputmax - inputmin outputrange = outputmax - outputmin OutputMatrix = (inputmatrix-inputmin) * outputrange/inputrange + outputmin return OutputMatrix + + +def simplified_solis(apparent_elevation, aod700=0.1, precipitable_water=1., + pressure=101325., dni_extra=1364.): + """ + Calculate the clear sky GHI, DNI, and DHI according to the + simplified Solis model [1]_. + + Reference [1]_ describes the accuracy of the model as being 15, 20, + and 18 W/m^2 for the beam, global, and diffuse components. Reference + [2]_ provides comparisons with other clear sky models. + + Parameters + ---------- + apparent_elevation: numeric + The apparent elevation of the sun above the horizon (deg). + + aod700: numeric + The aerosol optical depth at 700 nm (unitless). + Algorithm derived for values between 0 and 0.45. + + precipitable_water: numeric + The precipitable water of the atmosphere (cm). + Algorithm derived for values between 0.2 and 10 cm. + Values less than 0.2 will be assumed to be equal to 0.2. + + pressure: numeric + The atmospheric pressure (Pascals). + Algorithm derived for altitudes between sea level and 7000 m, + or 101325 and 41000 Pascals. + + dni_extra: numeric + Extraterrestrial irradiance. + + Returns + -------- + clearsky : DataFrame (if Series input) or OrderedDict of arrays + DataFrame/OrderedDict contains the columns/keys + ``'dhi', 'dni', 'ghi'``. + + The units of ``dni_extra`` determine the units of the output. + + References + ---------- + .. [1] P. Ineichen, "A broadband simplified version of the + Solis clear sky model," Solar Energy, 82, 758-762 (2008). + + .. [2] P. Ineichen, "Validation of models that estimate the clear + sky global and beam solar irradiance," Solar Energy, 132, + 332-344 (2016). + """ + + p = pressure + + w = precipitable_water + + # algorithm fails for pw < 0.2 + if np.isscalar(w): + w = 0.2 if w < 0.2 else w + else: + w = w.copy() + w[w < 0.2] = 0.2 + + # this algorithm is reasonably fast already, but it could be made + # faster by precalculating the powers of aod700, the log(p/p0), and + # the log(w) instead of repeating the calculations as needed in each + # function + + i0p = _calc_i0p(dni_extra, w, aod700, p) + + taub = _calc_taub(w, aod700, p) + b = _calc_b(w, aod700) + + taug = _calc_taug(w, aod700, p) + g = _calc_g(w, aod700) + + taud = _calc_taud(w, aod700, p) + d = _calc_d(w, aod700, p) + + # this prevents the creation of nans at night instead of 0s + # it's also friendly to scalar and series inputs + sin_elev = np.maximum(1.e-30, np.sin(np.radians(apparent_elevation))) + + dni = i0p * np.exp(-taub/sin_elev**b) + ghi = i0p * np.exp(-taug/sin_elev**g) * sin_elev + dhi = i0p * np.exp(-taud/sin_elev**d) + + irrads = OrderedDict() + irrads['ghi'] = ghi + irrads['dni'] = dni + irrads['dhi'] = dhi + + if isinstance(dni, pd.Series): + irrads = pd.DataFrame.from_dict(irrads) + + return irrads + + +def _calc_i0p(i0, w, aod700, p): + """Calculate the "enhanced extraterrestrial irradiance".""" + p0 = 101325. + io0 = 1.08 * w**0.0051 + i01 = 0.97 * w**0.032 + i02 = 0.12 * w**0.56 + i0p = i0 * (i02*aod700**2 + i01*aod700 + io0 + 0.071*np.log(p/p0)) + + return i0p + + +def _calc_taub(w, aod700, p): + """Calculate the taub coefficient""" + p0 = 101325. + tb1 = 1.82 + 0.056*np.log(w) + 0.0071*np.log(w)**2 + tb0 = 0.33 + 0.045*np.log(w) + 0.0096*np.log(w)**2 + tbp = 0.0089*w + 0.13 + + taub = tb1*aod700 + tb0 + tbp*np.log(p/p0) + + return taub + + +def _calc_b(w, aod700): + """Calculate the b coefficient.""" + + b1 = 0.00925*aod700**2 + 0.0148*aod700 - 0.0172 + b0 = -0.7565*aod700**2 + 0.5057*aod700 + 0.4557 + + b = b1 * np.log(w) + b0 + + return b + + +def _calc_taug(w, aod700, p): + """Calculate the taug coefficient""" + p0 = 101325. + tg1 = 1.24 + 0.047*np.log(w) + 0.0061*np.log(w)**2 + tg0 = 0.27 + 0.043*np.log(w) + 0.0090*np.log(w)**2 + tgp = 0.0079*w + 0.1 + taug = tg1*aod700 + tg0 + tgp*np.log(p/p0) + + return taug + + +def _calc_g(w, aod700): + """Calculate the g coefficient.""" + + g = -0.0147*np.log(w) - 0.3079*aod700**2 + 0.2846*aod700 + 0.3798 + + return g + + +def _calc_taud(w, aod700, p): + """Calculate the taud coefficient.""" + + # isscalar tests needed to ensure that the arrays will have the + # right shape in the tds calculation. + # there's probably a better way to do this. + + if np.isscalar(w) and np.isscalar(aod700): + w = np.array([w]) + aod700 = np.array([aod700]) + elif np.isscalar(w): + w = np.full_like(aod700, w) + elif np.isscalar(aod700): + aod700 = np.full_like(w, aod700) + + aod700_mask = aod700 < 0.05 + aod700_mask = np.array([aod700_mask, ~aod700_mask], dtype=np.int) + + # create tuples of coefficients for + # aod700 < 0.05, aod700 >= 0.05 + td4 = 86*w - 13800, -0.21*w + 11.6 + td3 = -3.11*w + 79.4, 0.27*w - 20.7 + td2 = -0.23*w + 74.8, -0.134*w + 15.5 + td1 = 0.092*w - 8.86, 0.0554*w - 5.71 + td0 = 0.0042*w + 3.12, 0.0057*w + 2.94 + tdp = -0.83*(1+aod700)**(-17.2), -0.71*(1+aod700)**(-15.0) + + tds = (np.array([td0, td1, td2, td3, td4, tdp]) * aod700_mask).sum(axis=1) + + p0 = 101325. + taud = (tds[4]*aod700**4 + tds[3]*aod700**3 + tds[2]*aod700**2 + + tds[1]*aod700 + tds[0] + tds[5]*np.log(p/p0)) + + # be polite about matching the output type to the input type(s) + if len(taud) == 1: + taud = taud[0] + + return taud + + +def _calc_d(w, aod700, p): + """Calculate the d coefficient.""" + + p0 = 101325. + dp = 1/(18 + 152*aod700) + d = -0.337*aod700**2 + 0.63*aod700 + 0.116 + dp*np.log(p/p0) + + return d diff --git a/pvlib/location.py b/pvlib/location.py index f926244a58..9a69edc2ee 100644 --- a/pvlib/location.py +++ b/pvlib/location.py @@ -9,62 +9,60 @@ import pandas as pd import pytz -from pvlib import solarposition -from pvlib import clearsky -from pvlib import atmosphere +from pvlib import solarposition, clearsky, atmosphere, irradiance class Location(object): """ Location objects are convenient containers for latitude, longitude, - timezone, and altitude data associated with a particular + timezone, and altitude data associated with a particular geographic location. You can also assign a name to a location object. - - Location objects have two timezone attributes: - + + Location objects have two timezone attributes: + * ``tz`` is a IANA timezone string. * ``pytz`` is a pytz timezone object. - + Location objects support the print method. - + Parameters ---------- latitude : float. Positive is north of the equator. Use decimal degrees notation. - - longitude : float. + + longitude : float. Positive is east of the prime meridian. Use decimal degrees notation. - - tz : str, int, float, or pytz.timezone. - See + + tz : str, int, float, or pytz.timezone. + See http://en.wikipedia.org/wiki/List_of_tz_database_time_zones for a list of valid time zones. pytz.timezone objects will be converted to strings. ints and floats must be in hours from UTC. - - alitude : float. + + alitude : float. Altitude from sea level in meters. - - name : None or string. + + name : None or string. Sets the name attribute of the Location object. - + **kwargs Arbitrary keyword arguments. Included for compatibility, but not used. - + See also -------- pvsystem.PVSystem """ - + def __init__(self, latitude, longitude, tz='UTC', altitude=0, name=None, **kwargs): - + self.latitude = latitude self.longitude = longitude - + if isinstance(tz, str): self.tz = tz self.pytz = pytz.timezone(tz) @@ -76,36 +74,34 @@ def __init__(self, latitude, longitude, tz='UTC', altitude=0, self.pytz = pytz.FixedOffset(tz*60) else: raise TypeError('Invalid tz specification') - + self.altitude = altitude - + self.name = name - + # needed for tying together Location and PVSystem in LocalizedPVSystem # if LocalizedPVSystem signature is reversed # super(Location, self).__init__(**kwargs) - - - + def __repr__(self): return ('{}: latitude={}, longitude={}, tz={}, altitude={}' - .format(self.name, self.latitude, self.longitude, + .format(self.name, self.latitude, self.longitude, self.tz, self.altitude)) - - + + @classmethod def from_tmy(cls, tmy_metadata, tmy_data=None, **kwargs): """ - Create an object based on a metadata + Create an object based on a metadata dictionary from tmy2 or tmy3 data readers. - + Parameters ---------- tmy_metadata : dict Returned from tmy.readtmy2 or tmy.readtmy3 tmy_data : None or DataFrame Optionally attach the TMY data to this object. - + Returns ------- Location object (or the child class of Location that you @@ -113,28 +109,28 @@ def from_tmy(cls, tmy_metadata, tmy_data=None, **kwargs): """ # not complete, but hopefully you get the idea. # might need code to handle the difference between tmy2 and tmy3 - + # determine if we're dealing with TMY2 or TMY3 data tmy2 = tmy_metadata.get('City', False) - + latitude = tmy_metadata['latitude'] longitude = tmy_metadata['longitude'] - + if tmy2: name = tmy_metadata['City'] - else: + else: name = tmy_metadata['Name'] - + tz = tmy_metadata['TZ'] altitude = tmy_metadata['altitude'] new_object = cls(latitude, longitude, tz=tz, altitude=altitude, name=name, **kwargs) - + # not sure if this should be assigned regardless of input. if tmy_data is not None: new_object.tmy_data = tmy_data - + return new_object @@ -143,7 +139,7 @@ def get_solarposition(self, times, pressure=None, temperature=12, """ Uses the :py:func:`solarposition.get_solarposition` function to calculate the solar zenith, azimuth, etc. at this location. - + Parameters ---------- times : DatetimeIndex @@ -153,7 +149,7 @@ def get_solarposition(self, times, pressure=None, temperature=12, temperature : None, float, or array-like kwargs passed to :py:func:`solarposition.get_solarposition` - + Returns ------- solar_position : DataFrame @@ -175,22 +171,24 @@ def get_clearsky(self, times, model='ineichen', **kwargs): """ Calculate the clear sky estimates of GHI, DNI, and/or DHI at this location. - + Parameters ---------- times : DatetimeIndex - + model : str - The clear sky model to use. - - kwargs passed to the relevant function(s). - + The clear sky model to use. Must be one of + 'ineichen', 'haurwitz', 'simplified_solis'. + + kwargs passed to the relevant functions. Climatological values + are assumed in many cases. See code for details. + Returns ------- clearsky : DataFrame Column names are: ``ghi, dni, dhi``. """ - + if model == 'ineichen': cs = clearsky.ineichen(times, latitude=self.latitude, longitude=self.longitude, @@ -199,18 +197,42 @@ def get_clearsky(self, times, model='ineichen', **kwargs): elif model == 'haurwitz': solpos = self.get_solarposition(times, **kwargs) cs = clearsky.haurwitz(solpos['apparent_zenith']) + elif model == 'simplified_solis': + + # these try/excepts define default values that are only + # evaluated if necessary. ineichen does some of this internally + try: + dni_extra = kwargs.pop('dni_extra') + except KeyError: + dni_extra = irradiance.extraradiation(times.dayofyear) + + try: + pressure = kwargs.pop('pressure') + except KeyError: + pressure = atmosphere.alt2pres(self.altitude) + + try: + apparent_elevation = kwargs.pop('apparent_elevation') + except KeyError: + solpos = self.get_solarposition( + times, pressure=pressure, **kwargs) + apparent_elevation = solpos['apparent_elevation'] + + cs = clearsky.simplified_solis( + apparent_elevation, pressure=pressure, dni_extra=dni_extra, + **kwargs) else: raise ValueError('{} is not a valid clear sky model' .format(model)) return cs - - + + def get_airmass(self, times=None, solar_position=None, model='kastenyoung1989'): """ Calculate the relative and absolute airmass. - + Automatically chooses zenith or apparant zenith depending on the selected model. @@ -222,7 +244,7 @@ def get_airmass(self, times=None, solar_position=None, DataFrame with with columns 'apparent_zenith', 'zenith'. model : str Relative airmass model - + Returns ------- airmass : DataFrame @@ -250,4 +272,3 @@ def get_airmass(self, times=None, solar_position=None, airmass['airmass_absolute'] = airmass_absolute return airmass - \ No newline at end of file diff --git a/pvlib/solarposition.py b/pvlib/solarposition.py index cce1ca66e5..188132166d 100644 --- a/pvlib/solarposition.py +++ b/pvlib/solarposition.py @@ -239,7 +239,7 @@ def _spa_python_import(how): def spa_python(time, latitude, longitude, altitude=0, pressure=101325, temperature=12, delta_t=None, - atmos_refract=None, how='numpy', numthreads=4): + atmos_refract=None, how='numpy', numthreads=4, **kwargs): """ Calculate the solar position using a python implementation of the NREL SPA algorithm described in [1]. diff --git a/pvlib/test/test_clearsky.py b/pvlib/test/test_clearsky.py index b6b50876d4..45e1957dc8 100644 --- a/pvlib/test/test_clearsky.py +++ b/pvlib/test/test_clearsky.py @@ -1,11 +1,13 @@ import logging pvl_logger = logging.getLogger('pvlib') +from collections import OrderedDict + import numpy as np import pandas as pd from nose.tools import raises -from numpy.testing import assert_almost_equal +from numpy.testing import assert_almost_equal, assert_allclose from pandas.util.testing import assert_frame_equal, assert_series_equal from pvlib.location import Location @@ -148,3 +150,217 @@ def test_haurwitz(): columns=['ghi'], index=times_localized) out = clearsky.haurwitz(ephem_data['zenith']) assert_frame_equal(expected, out) + + +def test_simplified_solis_series_elevation(): + expected = pd.DataFrame( + np.array([[ 0. , 0. , 0. ], + [ 0. , 0. , 0. ], + [ 377.80060035, 79.91931339, 42.77453223], + [ 869.47538184, 706.37903999, 110.05635962], + [ 958.89448856, 1062.44821373, 129.02349236], + [ 913.3209839 , 860.48978599, 118.94598678], + [ 634.01066762, 256.00505836, 72.18396705], + [ 0. , 0. , 0. ], + [ 0. , 0. , 0. ]]), + columns=['dni', 'ghi', 'dhi'], + index=times_localized) + expected = expected[['dhi', 'dni', 'ghi']] + + out = clearsky.simplified_solis(ephem_data['apparent_elevation']) + assert_frame_equal(expected, out) + + +def test_simplified_solis_scalar_elevation(): + expected = OrderedDict() + expected['ghi'] = 1064.653145 + expected['dni'] = 959.335463 + expected['dhi'] = 129.125602 + + out = clearsky.simplified_solis(80) + for k, v in expected.items(): + yield assert_allclose, expected[k], out[k] + + +def test_simplified_solis_scalar_neg_elevation(): + expected = OrderedDict() + expected['ghi'] = 0 + expected['dni'] = 0 + expected['dhi'] = 0 + + out = clearsky.simplified_solis(-10) + for k, v in expected.items(): + yield assert_allclose, expected[k], out[k] + + +def test_simplified_solis_series_elevation(): + expected = pd.DataFrame(np.array([[959.335463, 1064.653145, 129.125602]]), + columns=['dni', 'ghi', 'dhi']) + expected = expected[['ghi', 'dni', 'dhi']] + + out = clearsky.simplified_solis(pd.Series(80)) + assert_frame_equal(expected, out) + + +def test_simplified_solis_dni_extra(): + expected = pd.DataFrame(np.array([[963.555414, 1069.33637, 129.693603]]), + columns=['dni', 'ghi', 'dhi']) + expected = expected[['ghi', 'dni', 'dhi']] + + out = clearsky.simplified_solis(80, dni_extra=pd.Series(1370)) + assert_frame_equal(expected, out) + + +def test_simplified_solis_pressure(): + expected = pd.DataFrame(np. + array([[ 964.26930718, 1067.96543669, 127.22841797], + [ 961.88811874, 1066.36847963, 128.1402539 ], + [ 959.58112234, 1064.81837558, 129.0304193 ]]), + columns=['dni', 'ghi', 'dhi']) + expected = expected[['ghi', 'dni', 'dhi']] + + out = clearsky.simplified_solis( + 80, pressure=pd.Series([95000, 98000, 101000])) + assert_frame_equal(expected, out) + + +def test_simplified_solis_aod700(): + expected = pd.DataFrame(np. + array([[ 1056.61710493, 1105.7229086 , 64.41747323], + [ 1007.50558875, 1085.74139063, 102.96233698], + [ 959.3354628 , 1064.65314509, 129.12560167], + [ 342.45810926, 638.63409683, 77.71786575], + [ 55.24140911, 7.5413313 , 0. ]]), + columns=['dni', 'ghi', 'dhi']) + expected = expected[['ghi', 'dni', 'dhi']] + + aod700 = pd.Series([0.0, 0.05, 0.1, 1, 10]) + out = clearsky.simplified_solis(80, aod700=aod700) + assert_frame_equal(expected, out) + + +def test_simplified_solis_precipitable_water(): + expected = pd.DataFrame(np. + array([[ 1001.15353307, 1107.84678941, 128.58887606], + [ 1001.15353307, 1107.84678941, 128.58887606], + [ 983.51027357, 1089.62306672, 129.08755996], + [ 959.3354628 , 1064.65314509, 129.12560167], + [ 872.02335029, 974.18046717, 125.63581346]]), + columns=['dni', 'ghi', 'dhi']) + expected = expected[['ghi', 'dni', 'dhi']] + + out = clearsky.simplified_solis( + 80, precipitable_water=pd.Series([0.0, 0.2, 0.5, 1.0, 5.0])) + assert_frame_equal(expected, out) + + +def test_simplified_solis_small_scalar_pw(): + + expected = OrderedDict() + expected['ghi'] = 1107.84678941 + expected['dni'] = 1001.15353307 + expected['dhi'] = 128.58887606 + + out = clearsky.simplified_solis(80, precipitable_water=0.1) + for k, v in expected.items(): + yield assert_allclose, expected[k], out[k] + + +def test_simplified_solis_return_arrays(): + expected = OrderedDict() + + expected['ghi'] = np.array([[ 1148.40081325, 913.42330823], + [ 965.48550828, 760.04527609]]) + + expected['dni'] = np.array([[ 1099.25706525, 656.24601381], + [ 915.31689149, 530.31697378]]) + + expected['dhi'] = np.array([[ 64.1063074 , 254.6186615 ], + [ 62.75642216, 232.21931597]]) + + aod700 = np.linspace(0, 0.5, 2) + precipitable_water = np.linspace(0, 10, 2) + + aod700, precipitable_water = np.meshgrid(aod700, precipitable_water) + + out = clearsky.simplified_solis(80, aod700, precipitable_water) + + for k, v in expected.items(): + yield assert_allclose, expected[k], out[k] + + +def test_simplified_solis_nans_arrays(): + + # construct input arrays that each have 1 nan offset from each other, + # the last point is valid for all arrays + + length = 6 + + apparent_elevation = np.full(length, 80.) + apparent_elevation[0] = np.nan + + aod700 = np.full(length, 0.1) + aod700[1] = np.nan + + precipitable_water = np.full(length, 0.5) + precipitable_water[2] = np.nan + + pressure = np.full(length, 98000.) + pressure[3] = np.nan + + dni_extra = np.full(length, 1370.) + dni_extra[4] = np.nan + + expected = OrderedDict() + expected['ghi'] = np.full(length, np.nan) + expected['dni'] = np.full(length, np.nan) + expected['dhi'] = np.full(length, np.nan) + + expected['ghi'][length-1] = 1096.022736 + expected['dni'][length-1] = 990.306854 + expected['dhi'][length-1] = 128.664594 + + out = clearsky.simplified_solis(apparent_elevation, aod700, + precipitable_water, pressure, dni_extra) + + for k, v in expected.items(): + yield assert_allclose, expected[k], out[k] + + +def test_simplified_solis_nans_series(): + + # construct input arrays that each have 1 nan offset from each other, + # the last point is valid for all arrays + + length = 6 + + apparent_elevation = pd.Series(np.full(length, 80.)) + apparent_elevation[0] = np.nan + + aod700 = np.full(length, 0.1) + aod700[1] = np.nan + + precipitable_water = np.full(length, 0.5) + precipitable_water[2] = np.nan + + pressure = np.full(length, 98000.) + pressure[3] = np.nan + + dni_extra = np.full(length, 1370.) + dni_extra[4] = np.nan + + expected = OrderedDict() + expected['ghi'] = np.full(length, np.nan) + expected['dni'] = np.full(length, np.nan) + expected['dhi'] = np.full(length, np.nan) + + expected['ghi'][length-1] = 1096.022736 + expected['dni'][length-1] = 990.306854 + expected['dhi'][length-1] = 128.664594 + + expected = pd.DataFrame.from_dict(expected) + + out = clearsky.simplified_solis(apparent_elevation, aod700, + precipitable_water, pressure, dni_extra) + + assert_frame_equal(expected, out) diff --git a/pvlib/test/test_location.py b/pvlib/test/test_location.py index b2a9163d7b..42deab3038 100644 --- a/pvlib/test/test_location.py +++ b/pvlib/test/test_location.py @@ -79,6 +79,101 @@ def test_get_clearsky_haurwitz(): assert_frame_equal(expected, clearsky) +def test_get_clearsky_simplified_solis(): + tus = Location(32.2, -111, 'US/Arizona', 700, 'Tucson') + times = pd.DatetimeIndex(start='20160101T0600-0700', + end='20160101T1800-0700', + freq='3H') + clearsky = tus.get_clearsky(times, model='simplified_solis') + expected = pd.DataFrame(data=np. + array([[ 0. , 0. , 0. ], + [ 70.00146271, 638.01145669, 236.71136245], + [ 101.69729217, 852.51950946, 577.1117803 ], + [ 86.1679965 , 755.98048017, 385.59586091], + [ 0. , 0. , 0. ]]), + columns=['dhi', 'dni', 'ghi'], + index=times) + expected = expected[['ghi', 'dni', 'dhi']] + assert_frame_equal(expected, clearsky) + + +def test_get_clearsky_simplified_solis_apparent_elevation(): + tus = Location(32.2, -111, 'US/Arizona', 700, 'Tucson') + times = pd.DatetimeIndex(start='20160101T0600-0700', + end='20160101T1800-0700', + freq='3H') + apparent_elevation = pd.Series(80, index=times) + clearsky = tus.get_clearsky(times, model='simplified_solis', + apparent_elevation=apparent_elevation) + expected = pd.DataFrame(data=np. + array([[ 131.3124497 , 1001.14754036, 1108.14147919], + [ 131.3124497 , 1001.14754036, 1108.14147919], + [ 131.3124497 , 1001.14754036, 1108.14147919], + [ 131.3124497 , 1001.14754036, 1108.14147919], + [ 131.3124497 , 1001.14754036, 1108.14147919]]), + columns=['dhi', 'dni', 'ghi'], + index=times) + expected = expected[['ghi', 'dni', 'dhi']] + assert_frame_equal(expected, clearsky) + + +def test_get_clearsky_simplified_solis_dni_extra(): + tus = Location(32.2, -111, 'US/Arizona', 700, 'Tucson') + times = pd.DatetimeIndex(start='20160101T0600-0700', + end='20160101T1800-0700', + freq='3H') + clearsky = tus.get_clearsky(times, model='simplified_solis', + dni_extra=1370) + expected = pd.DataFrame(data=np. + array([[ 0. , 0. , 0. ], + [ 67.82281485, 618.15469596, 229.34422063], + [ 98.53217848, 825.98663808, 559.15039353], + [ 83.48619937, 732.45218243, 373.59500313], + [ 0. , 0. , 0. ]]), + columns=['dhi', 'dni', 'ghi'], + index=times) + expected = expected[['ghi', 'dni', 'dhi']] + assert_frame_equal(expected, clearsky) + + +def test_get_clearsky_simplified_solis_pressure(): + tus = Location(32.2, -111, 'US/Arizona', 700, 'Tucson') + times = pd.DatetimeIndex(start='20160101T0600-0700', + end='20160101T1800-0700', + freq='3H') + clearsky = tus.get_clearsky(times, model='simplified_solis', + pressure=95000) + expected = pd.DataFrame(data=np. + array([[ 0. , 0. , 0. ], + [ 70.20556637, 635.53091983, 236.17716435], + [ 102.08954904, 850.49502085, 576.28465815], + [ 86.46561686, 753.70744638, 384.90537859], + [ 0. , 0. , 0. ]]), + columns=['dhi', 'dni', 'ghi'], + index=times) + expected = expected[['ghi', 'dni', 'dhi']] + assert_frame_equal(expected, clearsky) + + +def test_get_clearsky_simplified_solis_aod_pw(): + tus = Location(32.2, -111, 'US/Arizona', 700, 'Tucson') + times = pd.DatetimeIndex(start='20160101T0600-0700', + end='20160101T1800-0700', + freq='3H') + clearsky = tus.get_clearsky(times, model='simplified_solis', + aod700=0.25, precipitable_water=2.) + expected = pd.DataFrame(data=np. + array([[ 0. , 0. , 0. ], + [ 85.77821205, 374.58084365, 179.48483117], + [ 143.52743364, 625.91745295, 490.06254157], + [ 114.63275842, 506.52275195, 312.24711495], + [ 0. , 0. , 0. ]]), + columns=['dhi', 'dni', 'ghi'], + index=times) + expected = expected[['ghi', 'dni', 'dhi']] + assert_frame_equal(expected, clearsky) + + @raises(ValueError) def test_get_clearsky_valueerror(): tus = Location(32.2, -111, 'US/Arizona', 700, 'Tucson') @@ -160,10 +255,10 @@ def test_get_airmass_valueerror(): end='20160101T1800-0700', freq='3H') clearsky = tus.get_airmass(times, model='invalid_model') - + def test_Location___repr__(): tus = Location(32.2, -111, 'US/Arizona', 700, 'Tucson') assert tus.__repr__()==('Tucson: latitude=32.2, longitude=-111, '+ 'tz=US/Arizona, altitude=700') - - + +