From 33830cb15560b42eb8fde29619273dd712de36a2 Mon Sep 17 00:00:00 2001 From: "Marc A. Anoma" Date: Wed, 12 Sep 2018 15:16:44 +0900 Subject: [PATCH] Simplify timeseries calculation API (#14) * Remove "save_segments" arguments as not needed anymore * Compact formatting in timeseries calculation function * Fix issue with test_serial_perez that wasn't caught in previous PR * Improving test for perez luminance calculation * Major update of timeseries calculation API * all tests are now passing * now passing individual arguments instead of dataframes with all arguments * broke custom perez calculation from pvarray timeseries calculation: this will allow the users to use one or the other separately * next step will be to remove the "simple" calculation mode * Add pytest-mock to circleci config for testing * Removed the specific functions related to "isotropic" calcs * now only using the perez functions, even when using isotropic approach (just setting circumsolar and horizon to 0) * Reorganizing package: seperate timeseries fns from plotting ones * Update docstrings in timeseries and plot modules --- .circleci/config.yml | 2 +- pvfactors/__init__.py | 28 + pvfactors/plot.py | 117 ++++ pvfactors/pvarray.py | 33 - pvfactors/tests/conftest.py | 8 + pvfactors/tests/test_array_geometry.py | 2 +- pvfactors/tests/test_circumsolar_shading.py | 14 +- pvfactors/tests/test_discretization.py | 10 +- .../file_test_df_perez_luminance.csv | 5 + .../test_files/file_test_serial_perez.csv | 2 +- pvfactors/tests/test_multiprocessing.py | 12 +- .../tests/test_serial_perez_calculation.py | 21 +- pvfactors/tests/test_tools.py | 126 +++- pvfactors/tests/test_view_factor.py | 19 +- pvfactors/timeseries.py | 487 ++++++++++++++ pvfactors/tools.py | 603 ------------------ setup.py | 2 +- 17 files changed, 797 insertions(+), 694 deletions(-) create mode 100644 pvfactors/plot.py create mode 100644 pvfactors/tests/test_files/file_test_df_perez_luminance.csv create mode 100644 pvfactors/timeseries.py delete mode 100644 pvfactors/tools.py diff --git a/.circleci/config.yml b/.circleci/config.yml index 43b7b872..cf9ab5e2 100644 --- a/.circleci/config.yml +++ b/.circleci/config.yml @@ -23,7 +23,7 @@ jobs: - v1-dep-{{ .Branch }}- - v1-dep-master- - v1-dep- - - run: sudo pip install pytest==3.2.2 numpy==1.13.3 scipy==0.19.1 pandas==0.23.3 shapely==1.6.1 pvlib==0.5.0 future==0.16.0 six==1.11.0 + - run: sudo pip install pytest==3.2.2 numpy==1.13.3 scipy==0.19.1 pandas==0.23.3 shapely==1.6.1 pvlib==0.5.0 future==0.16.0 six==1.11.0 pytest-mock==1.10.0 - save_cache: key: v1-dep-{{ .Branch }}-{{ epoch }} paths: diff --git a/pvfactors/__init__.py b/pvfactors/__init__.py index 3feddfeb..294ab09e 100644 --- a/pvfactors/__init__.py +++ b/pvfactors/__init__.py @@ -1,6 +1,7 @@ # -*- coding: utf-8 -*- import logging +import sys logging.basicConfig() from pvfactors.version import __version__ @@ -16,3 +17,30 @@ class PVFactorsEdgePointDoesNotExist(Exception): class PVFactorsArrayUpdateException(Exception): pass + + +# Define function used for progress bar when running long simulations +# Borrowed from: https://gist.github.com/aubricus +def print_progress(iteration, total, prefix='', suffix='', decimals=1, + bar_length=100): + """ + Call in a loop to create terminal progress bar + @params: + iteration - Required : current iteration (Int) + total - Required : total iterations (Int) + prefix - Optional : prefix string (Str) + suffix - Optional : suffix string (Str) + decimals - Optional : positive number of decimals in percent + complete (Int) + bar_length - Optional : character length of bar (Int) + """ + format_str = "{0:." + str(decimals) + "f}" + percents = format_str.format(100 * (iteration / float(total))) + filled_length = int(round(bar_length * iteration / float(total))) + bar = '█' * filled_length + '-' * (bar_length - filled_length) + sys.stdout.write('\r%s |%s| %s%s %s' % (prefix, bar, percents, '%', + suffix)), + sys.stdout.flush() + if iteration == total: + sys.stdout.write('\n') + sys.stdout.flush() diff --git a/pvfactors/plot.py b/pvfactors/plot.py new file mode 100644 index 00000000..0ea71335 --- /dev/null +++ b/pvfactors/plot.py @@ -0,0 +1,117 @@ +# -*- coding: utf-8 -*- + +import numpy as np +from pvfactors import logging + +LOGGER = logging.getLogger(__name__) +LOGGER.setLevel(logging.INFO) + +# Define colors used for plotting the 2D arrays +COLOR_dic = { + 'i': '#FFBB33', + 's': '#A7A49D', + 't': '#6699cc', + 'pvrow_illum': '#6699cc', + 'pvrow_shaded': '#ff0000', + 'ground_shaded': '#A7A49D', + 'ground_illum': '#FFBB33' +} + + +def plot_array_from_registry(ax, registry, line_types_selected=None, + fontsize=20): + """ + Plot a 2D PV array using the ``shapely`` geometry objects located in + a :class:`pvarray.Array` surface registry. + + :param matplotlib.axes.Axes ax: axes to use for the plot + :param pd.DataFrame registry: registry containing geometries to plot + :param list line_types_selected: parameter used to select a subset of + 'line_type' to plot; e.g. 'pvrow' or 'ground' + :return: None (``ax`` is updated) + """ + + registry = registry.copy() + registry.loc[:, 'color'] = ( + registry.line_type.values + '_' + + np.where(registry.shaded.values, 'shaded', 'illum')) + # TODO: distance may not exist + if line_types_selected: + for line_type in line_types_selected: + surface_reg_selected = registry.loc[ + registry.line_type == line_type, :] + for index, row in surface_reg_selected.iterrows(): + LOGGER.debug("Plotting %s", row['line_type']) + plot_coords(ax, row['geometry']) + plot_bounds(ax, row['geometry']) + plot_line(ax, row['geometry'], row['style'], + row['shading_type']) + else: + for index, row in registry.iterrows(): + LOGGER.debug("Plotting %s", row['line_type']) + plot_coords(ax, row['geometry']) + plot_bounds(ax, row['geometry']) + plot_line(ax, row['geometry'], row['style'], row['color']) + + ax.axis('equal') + ax.set_xlabel("x [m]", fontsize=fontsize) + ax.set_ylabel("y [m]", fontsize=fontsize) + + +def plot_pvarray(ax, pvarray, line_types_selected=None, fontsize=20): + """ + Plot a 2D PV array from a :class:`pvarray.Array` using its + :attr:`pvarray.Array.surface_registry`. + + :param ax: :class:`matplotlib.axes.Axes` object to use for the plot + :param pvarray: object containing the surface registry as attribute + :type pvarray: :class:`pvarray.Array` + :param list line_types_selected: parameter used to select a subset of + 'line_type' to plot; e.g. 'pvrow' or 'ground' + :return: None (``ax`` is updated) + """ + + # FIXME: repeating code from plot_line_registry + surface_registry = pvarray.surface_registry.copy() + plot_array_from_registry(ax, surface_registry, + line_types_selected=line_types_selected) + + # Plot details + distance = pvarray.pvrow_distance + height = pvarray.pvrow_height + n_pvrows = pvarray.n_pvrows + ax.set_xlim(- 0.5 * distance, (n_pvrows - 0.5) * distance) + ax.set_ylim(-height, 2 * height) + ax.set_title("PV Array", fontsize=fontsize) + + +# Base functions used to plot the 2D array +def plot_coords(ax, ob): + try: + x, y = ob.xy + ax.plot(x, y, 'o', color='#999999', zorder=1) + except NotImplementedError: + for line in ob: + x, y = line.xy + ax.plot(x, y, 'o', color='#999999', zorder=1) + + +def plot_bounds(ax, ob): + # Check if shadow reduces to one point (for very specific sun alignment) + if len(ob.boundary) == 0: + x, y = ob.coords[0] + else: + x, y = zip(*list((p.x, p.y) for p in ob.boundary)) + ax.plot(x, y, 'o', color='#000000', zorder=1) + + +def plot_line(ax, ob, line_style, line_color): + try: + x, y = ob.xy + ax.plot(x, y, color=COLOR_dic[line_color], ls=line_style, alpha=0.7, + linewidth=3, solid_capstyle='round', zorder=2) + except NotImplementedError: + for line in ob: + x, y = line.xy + ax.plot(x, y, color=COLOR_dic[line_color], ls=line_style, + alpha=0.7, linewidth=3, solid_capstyle='round', zorder=2) diff --git a/pvfactors/pvarray.py b/pvfactors/pvarray.py index bffc376a..112a9825 100644 --- a/pvfactors/pvarray.py +++ b/pvfactors/pvarray.py @@ -511,39 +511,6 @@ def update_reflectivity_matrix(self): self.inv_reflectivity_matrix = np.diag( list(1. / self.surface_registry.reflectivity.values) + [1]) - def calculate_radiosities_simple(self, solar_zenith, solar_azimuth, - array_tilt, array_azimuth, dni, dhi): - """ - Solve linear system of equations to calculate radiosity terms based on - the specified inputs and assuming that the whole sky dome is isotropic - (no diffuse sky dome decomposition) - - :param float solar_zenith: zenith angle of the sun [degrees] - :param float solar_azimuth: azimuth angle of the sun [degrees] - :param float array_tilt: tilt angle of the whole array. All PV rows must - have the same tilt angle [degrees] - :param float array_azimuth: azimuth angle of the whole array. All PV - rows must have the same azimuth angle [degrees] - :param float dni: direct normal irradiance [W/m2] - :param float dhi: diffuse horizontal irradiance [W/m2] - :return: None; updating :attr:`surface_registry` - """ - # Update the array configuration - self.update_view_factors(solar_zenith, solar_azimuth, array_tilt, - array_azimuth) - self.update_irradiance_terms_simple(solar_zenith, solar_azimuth, - array_tilt, array_azimuth, dni, dhi) - self.update_reflectivity_matrix() - - # Do calculation - a_mat = self.inv_reflectivity_matrix - self.vf_matrix - q0 = linalg.solve(a_mat, self.irradiance_terms) - qinc = np.dot(self.vf_matrix, q0) + self.irradiance_terms - # Assign to surfaces - self.surface_registry.loc[:, 'q0'] = q0[:-1] - self.surface_registry.loc[:, 'qinc'] = qinc[:-1] - self.calculate_sky_and_reflection_components() - def calculate_radiosities_perez( self, solar_zenith, solar_azimuth, array_tilt, array_azimuth, dni, luminance_isotropic, luminance_circumsolar, poa_horizon, diff --git a/pvfactors/tests/conftest.py b/pvfactors/tests/conftest.py index e0b516a1..019e7568 100644 --- a/pvfactors/tests/conftest.py +++ b/pvfactors/tests/conftest.py @@ -33,3 +33,11 @@ def df_segments(): df_segments = pd.read_csv(fp, header=[0, 1], index_col=0) df_segments.index = pd.to_datetime(df_segments.index) yield df_segments + + +@pytest.fixture(scope='function') +def df_perez_luminance(): + """ Example of df_segments to be used for tests """ + fp = os.path.join(DIR_TEST_DATA, 'file_test_df_perez_luminance.csv') + df_perez_luminance = pd.read_csv(fp, header=[0], index_col=0) + yield df_perez_luminance diff --git a/pvfactors/tests/test_array_geometry.py b/pvfactors/tests/test_array_geometry.py index e6cc3fac..e022e1fb 100644 --- a/pvfactors/tests/test_array_geometry.py +++ b/pvfactors/tests/test_array_geometry.py @@ -51,7 +51,7 @@ def test_plotting(): is_ci = os.environ.get('CI', False) if not is_ci: import matplotlib.pyplot as plt - from pvfactors.tools import plot_pvarray + from pvfactors.plot import plot_pvarray # Create array where sun vector is in the direction of the modules arguments = { 'n_pvrows': 3, diff --git a/pvfactors/tests/test_circumsolar_shading.py b/pvfactors/tests/test_circumsolar_shading.py index 189775c1..c3c76723 100644 --- a/pvfactors/tests/test_circumsolar_shading.py +++ b/pvfactors/tests/test_circumsolar_shading.py @@ -6,8 +6,8 @@ """ from pvfactors.pvcore import calculate_circumsolar_shading -from pvfactors.pvarray import Array -from pvfactors.tools import calculate_radiosities_serially_perez +from pvfactors.timeseries import (calculate_radiosities_serially_perez, + breakup_df_inputs) import pandas as pd import os import numpy as np @@ -61,17 +61,17 @@ def test_serial_circumsolar_shading_calculation(): 'calculate_front_circ_horizon_shading': True, 'circumsolar_model': 'gaussian' } - save = (1, 'front') # Load inputs for the serial calculation test_file = os.path.join( TEST_DATA, 'file_test_serial_circumsolar_shading_calculation.csv') df_inputs = pd.read_csv(test_file, index_col=0) df_inputs.index = pd.DatetimeIndex(df_inputs.index) - - # Create shapely PV array - array = Array(**arguments) + (timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, dni, dhi) = breakup_df_inputs(df_inputs) # Run the calculation for functional testing df_registries, df_inputs_perez = ( - calculate_radiosities_serially_perez((arguments, df_inputs, save)) + calculate_radiosities_serially_perez((arguments, timestamps, array_tilt, + array_azimuth, solar_zenith, + solar_azimuth, dni, dhi)) ) diff --git a/pvfactors/tests/test_discretization.py b/pvfactors/tests/test_discretization.py index 806f4964..95d85f7b 100644 --- a/pvfactors/tests/test_discretization.py +++ b/pvfactors/tests/test_discretization.py @@ -60,6 +60,10 @@ def test_consistent_qinc(): # Run a calculation for the given configuration dni = 1e3 dhi = 1e2 + luminance_isotropic = dhi + luminance_circumsolar = 0. + poa_horizon = 0. + poa_circumsolar = 0. solar_zenith = 20. solar_azimuth = 180. @@ -67,8 +71,10 @@ def test_consistent_qinc(): array_tilt = 20. array_azimuth = 180. - array.calculate_radiosities_simple(solar_zenith, solar_azimuth, array_tilt, - array_azimuth, dni, dhi) + array.calculate_radiosities_perez(solar_zenith, solar_azimuth, array_tilt, + array_azimuth, dni, luminance_isotropic, + luminance_circumsolar, poa_horizon, + poa_circumsolar) # Compare to expected values expected_qinc = np.array([ diff --git a/pvfactors/tests/test_files/file_test_df_perez_luminance.csv b/pvfactors/tests/test_files/file_test_df_perez_luminance.csv new file mode 100644 index 00000000..e245585e --- /dev/null +++ b/pvfactors/tests/test_files/file_test_df_perez_luminance.csv @@ -0,0 +1,5 @@ +datetime,solar_zenith,solar_azimuth,array_tilt,array_azimuth,dni,dhi,poa_isotropic,poa_circumsolar,poa_horizon,vf_horizon,vf_circumsolar,vf_isotropic,luminance_horizon,luminance_circumsolar,luminance_isotropic,poa_total_diffuse +2015-10-02 12:40:00.001200,42.54210405,174.3715553,5.142972856,90,1000,100,53.9049589579,46.1721986001,1.8417661189,0.0896413203484,1.00404216032,0.99798705648,20.5459503691,45.986314544,54.013685456,101.918923677 +2015-10-02 12:49:59.998800,42.41609807,178.0613184,1.77035709,90,1000,100,53.9393074992,46.0698063615,0.633882293345,0.0308936438071,1.00047755049,0.999761338734,20.5182107136,46.0478162043,53.9521837957,100.642996154 +2015-10-02 13:00:00.000000,42.41574779,181.7602662,1.607556458,270,1000,100,53.9413957759,46.0661176346,0.575604882308,0.0280534721322,1.00039373104,0.99980321195,20.518133356,46.0479871126,53.9520128874,100.583118293 +2015-10-02 13:10:00.001200,42.54106148,185.4501296,4.98118732,270,1000,100,53.9111795968,46.161162222,1.78395715039,0.0868286442933,1.00379105643,0.998111630694,20.545721575,45.9868235788,54.0131764212,101.856298969 diff --git a/pvfactors/tests/test_files/file_test_serial_perez.csv b/pvfactors/tests/test_files/file_test_serial_perez.csv index f3e5cd51..ac2ef32d 100644 --- a/pvfactors/tests/test_files/file_test_serial_perez.csv +++ b/pvfactors/tests/test_files/file_test_serial_perez.csv @@ -1,4 +1,4 @@ -timestamps,pvrow,side,term,value +timestamps,pvrow_index,surface_side,term,value 4/15/2017 10:26,0,back,circumsolar_term,0 4/15/2017 10:27,0,back,circumsolar_term,0 4/15/2017 10:26,0,front,circumsolar_term,51.02503788 diff --git a/pvfactors/tests/test_multiprocessing.py b/pvfactors/tests/test_multiprocessing.py index c6a6f652..1c18f1ef 100644 --- a/pvfactors/tests/test_multiprocessing.py +++ b/pvfactors/tests/test_multiprocessing.py @@ -4,7 +4,8 @@ Testing the multiprocessing calculation from the ``tools`` module """ -from pvfactors.tools import calculate_radiosities_parallel_perez +from pvfactors.timeseries import (calculate_radiosities_parallel_perez, + breakup_df_inputs) import os import pandas as pd @@ -43,7 +44,10 @@ def test_calculate_radiosities_parallel_perez(): df_inputs_simulation = df_inputs_simulation.iloc[:subset_idx, :] # Select number of processes n_processes = None + # Break up inputs + (timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, dni, dhi) = breakup_df_inputs(df_inputs_simulation) # Run calculation - results = calculate_radiosities_parallel_perez(arguments, - df_inputs_simulation, - n_processes=n_processes) + results = calculate_radiosities_parallel_perez( + arguments, timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, dni, dhi, n_processes=n_processes) diff --git a/pvfactors/tests/test_serial_perez_calculation.py b/pvfactors/tests/test_serial_perez_calculation.py index 2727fc61..1ee1702f 100644 --- a/pvfactors/tests/test_serial_perez_calculation.py +++ b/pvfactors/tests/test_serial_perez_calculation.py @@ -5,8 +5,9 @@ consistent results """ -from pvfactors.tools import (calculate_radiosities_serially_perez, - get_average_pvrow_outputs) +from pvfactors.timeseries import (calculate_radiosities_serially_perez, + get_average_pvrow_outputs, + breakup_df_inputs) import pandas as pd import numpy as np import os @@ -43,13 +44,18 @@ def test_serial_calculation(): idx_subset = 10 df_inputs_simulation = df_inputs_simulation.iloc[0:idx_subset, :] + # Break up inputs + (timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, dni, dhi) = breakup_df_inputs(df_inputs_simulation) + # Run calculation in 1 process only - (df_registries, _) = ( - calculate_radiosities_serially_perez((arguments, df_inputs_simulation)) - ) + df_registries, _ = calculate_radiosities_serially_perez( + (arguments, timestamps, + solar_zenith, solar_azimuth, + array_tilt, array_azimuth, dni, dhi)) # Format df_registries to get outputs - df_outputs = get_average_pvrow_outputs(df_registries) + df_outputs = get_average_pvrow_outputs(df_registries, include_shading=False) # Did the outputs remain consistent? test_results = values_are_consistent(df_outputs) @@ -74,7 +80,6 @@ def values_are_consistent(df_outputs): parse_dates=['timestamps']) expected_df_outputs['origin'] = 'expected' df_outputs = (df_outputs.assign(timestamps=lambda x: x.index) - .drop('shaded', axis=1) .reset_index(drop=True) .melt(id_vars=['timestamps'])) df_outputs['origin'] = 'calculated' @@ -89,7 +94,7 @@ def values_are_consistent(df_outputs): rtol = 1e-7 test_results = [] for name, group in grouped_comparison: - df_term = group.pivot_table(index=['timestamps', 'pvrow', 'side'], + df_term = group.pivot_table(index=['timestamps', 'pvrow_index', 'surface_side'], columns=['origin'], values='value') compare = (('calculated' in df_term.columns) & ('expected' in df_term.columns)) diff --git a/pvfactors/tests/test_tools.py b/pvfactors/tests/test_tools.py index 6e49d97c..0acdfc3b 100644 --- a/pvfactors/tests/test_tools.py +++ b/pvfactors/tests/test_tools.py @@ -3,15 +3,15 @@ """ Test some of the functions in the tools module """ - +import pytest from pvfactors.pvarray import Array -from pvfactors.tools import (calculate_radiosities_serially_simple, - perez_diffuse_luminance, - calculate_radiosities_serially_perez, - calculate_radiosities_parallel_perez, - get_average_pvrow_outputs, - get_bifacial_gain_outputs, - get_pvrow_segment_outputs) +from pvfactors.timeseries import (perez_diffuse_luminance, + calculate_radiosities_serially_perez, + calculate_radiosities_parallel_perez, + get_average_pvrow_outputs, + get_pvrow_segment_outputs, + breakup_df_inputs, + array_timeseries_calculate) import numpy as np import pandas as pd import os @@ -22,10 +22,15 @@ idx_slice = pd.IndexSlice -def test_calculate_radiosities_serially_simple(): +@pytest.fixture(scope='function') +def mock_array_timeseries_calculate(mocker): + return mocker.patch('pvfactors.timeseries.array_timeseries_calculate') + + +def test_array_calculate_timeseries(): """ - Check that the results of the radiosity calculation using the isotropic - diffuse sky stay consistent + Check that the timeseries results of the radiosity calculation using the isotropic + diffuse sky approach stay consistent """ # Simple sky and array configuration df_inputs = pd.DataFrame( @@ -33,8 +38,8 @@ def test_calculate_radiosities_serially_simple(): [[80., 0., 70., 180., 1e3, 1e2], [20., 180., 40., 180., 1e3, 1e2], [70.4407256, 248.08690811, 42.4337927, 270., 1000., 100.]]), - columns=['solar_zenith', 'solar_azimuth', 'array_tilt', 'array_azimuth', - 'dni', 'dhi'], + columns=['solar_zenith', 'solar_azimuth', 'array_tilt', + 'array_azimuth', 'dni', 'dhi'], index=[0, 1, 2] ) arguments = { @@ -43,11 +48,31 @@ def test_calculate_radiosities_serially_simple(): 'pvrow_width': 1., 'gcr': 0.3, } + array = Array(**arguments) - # Calculate irradiance terms - df_outputs, df_bifacial_gains = ( - calculate_radiosities_serially_simple(array, df_inputs)) + # Break up inputs + (timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, dni, dhi) = breakup_df_inputs(df_inputs) + + # Fill in the missing pieces + luminance_isotropic = dhi + luminance_circumsolar = np.zeros(len(timestamps)) + poa_horizon = np.zeros(len(timestamps)) + poa_circumsolar = np.zeros(len(timestamps)) + + # # Calculate irradiance terms + # df_outputs, df_bifacial_g = ( + # calculate_radiosities_serially_simple(array, df_inputs)) + + # Run timeseries calculation + df_registries = array_timeseries_calculate( + arguments, timestamps, solar_zenith, solar_azimuth, array_tilt, array_azimuth, + dni, luminance_isotropic, luminance_circumsolar, poa_horizon, poa_circumsolar) + + # Calculate surface averages for pvrows + df_outputs = get_average_pvrow_outputs(df_registries, values=['q0', 'qinc'], + include_shading=False) # Check that the outputs are as expected expected_outputs_array = np.array([ @@ -66,11 +91,36 @@ def test_calculate_radiosities_serially_simple(): [True, False, False]], dtype=object) tol = 1e-8 assert np.allclose(expected_outputs_array[:-1, :].astype(float), - df_outputs.values[:-1, :].astype(float), + df_outputs.values.T, atol=tol, rtol=0, equal_nan=True) -def test_perez_diffuse_luminance(): +def test_perez_diffuse_luminance(df_perez_luminance): + """ + Test that the calculation of luminance -- first step in using the vf model + with Perez -- is functional + """ + df_inputs_clearday = pd.read_csv(FILE_PATH) + df_inputs_clearday = df_inputs_clearday.set_index('datetime', drop=True) + df_inputs_clearday.index = (pd.DatetimeIndex(df_inputs_clearday.index) + .tz_localize('UTC').tz_convert('Etc/GMT+7') + .tz_localize(None)) + + # Break up inputs + (timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, dni, dhi) = breakup_df_inputs(df_inputs_clearday) + df_outputs = perez_diffuse_luminance(timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, dni, dhi) + + col_order = df_outputs.columns + tol = 1e-8 + assert np.allclose(df_outputs.values, + df_perez_luminance[col_order].values, + atol=0, rtol=tol) + + +def test_luminance_in_timeseries_calc(df_perez_luminance, + mock_array_timeseries_calculate): """ Test that the calculation of luminance -- first step in using the vf model with Perez -- is functional @@ -81,7 +131,19 @@ def test_perez_diffuse_luminance(): .tz_localize('UTC').tz_convert('Etc/GMT+7') .tz_localize(None)) - df_outputs = perez_diffuse_luminance(df_inputs_clearday) + # Break up inputs + (timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, dni, dhi) = breakup_df_inputs(df_inputs_clearday) + _, df_outputs = calculate_radiosities_serially_perez( + (None, timestamps, solar_zenith, solar_azimuth, + array_tilt, array_azimuth, + dni, dhi)) + + col_order = df_outputs.columns + tol = 1e-8 + assert np.allclose(df_outputs.values, + df_perez_luminance[col_order].values, + atol=0, rtol=tol) def test_save_all_outputs_calculate_perez(): @@ -118,19 +180,20 @@ def test_save_all_outputs_calculate_perez(): 'cut': [(1, 3, 'front')] } - # We want to save the results from the front side of tracker #2 (index 1) - save_segments = (1, 'front') - args = (arguments, df_inputs_clearday.iloc[:idx_subset], save_segments) + # Break up inputs + (timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, dni, dhi) = breakup_df_inputs( + df_inputs_clearday.iloc[:idx_subset]) + + args = (arguments, timestamps, solar_zenith, solar_azimuth, + array_tilt, array_azimuth, dni, dhi) # Run the serial calculation df_registries_serial, _ = ( calculate_radiosities_serially_perez(args)) df_registries_parallel, _ = ( - calculate_radiosities_parallel_perez( - arguments, df_inputs_clearday.iloc[:idx_subset], - save_segments=save_segments - )) + calculate_radiosities_parallel_perez(*args)) # Format the outputs df_outputs_segments_serial = get_pvrow_segment_outputs( @@ -159,16 +222,21 @@ def test_save_all_outputs_calculate_perez(): def test_get_average_pvrow_outputs(df_registries, df_outputs): """ Test that obtaining the correct format, using outputs from v011 """ - calc_df_outputs = get_average_pvrow_outputs(df_registries) + calc_df_outputs_num = get_average_pvrow_outputs(df_registries, + include_shading=False) + calc_df_outputs_shading = get_average_pvrow_outputs(df_registries, values=[], + include_shading=True) tol = 1e-8 # Compare numerical values assert np.allclose(df_outputs.iloc[:, 1:].values, - calc_df_outputs.iloc[:, 1:].values, + calc_df_outputs_num.values, atol=0, rtol=tol) + + array_is_shaded = calc_df_outputs_shading.sum(axis=1).values > 0 # Compare bool on shading assert np.array_equal(df_outputs.iloc[:, 0].values, - calc_df_outputs.iloc[:, 0].values) + array_is_shaded) def test_get_average_pvrow_segments(df_registries, df_segments): diff --git a/pvfactors/tests/test_view_factor.py b/pvfactors/tests/test_view_factor.py index 17ab5c2b..52e8c3f1 100644 --- a/pvfactors/tests/test_view_factor.py +++ b/pvfactors/tests/test_view_factor.py @@ -5,8 +5,9 @@ """ from pvfactors.pvarray import Array -from pvfactors.tools import (calculate_radiosities_serially_perez, - get_average_pvrow_outputs) +from pvfactors.timeseries import (calculate_radiosities_serially_perez, + get_average_pvrow_outputs, + breakup_df_inputs) import numpy as np import pandas as pd import os @@ -144,7 +145,12 @@ def test_negativevf_and_flatcasenoon(): df_inputs.index = pd.DatetimeIndex(df_inputs.index).tz_localize( 'UTC').tz_convert('US/Arizona') - args = (pvarray_parameters, df_inputs) + # Break up inputs + (timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, dni, dhi) = breakup_df_inputs(df_inputs) + + args = (pvarray_parameters, timestamps, solar_zenith, solar_azimuth, + array_tilt, array_azimuth, dni, dhi) df_registries, _ = calculate_radiosities_serially_perez(args) df_outputs = get_average_pvrow_outputs(df_registries) @@ -187,7 +193,12 @@ def test_back_surface_luminance(): df_inputs.index = pd.DatetimeIndex(df_inputs.index).tz_localize( 'UTC').tz_convert('US/Arizona') - args = (pvarray_parameters, df_inputs) + # Break up inputs + (timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, dni, dhi) = breakup_df_inputs(df_inputs) + + args = (pvarray_parameters, timestamps, solar_zenith, solar_azimuth, + array_tilt, array_azimuth, dni, dhi) df_registries, _ = calculate_radiosities_serially_perez(args) df_outputs = get_average_pvrow_outputs(df_registries) diff --git a/pvfactors/timeseries.py b/pvfactors/timeseries.py new file mode 100644 index 00000000..00ecec19 --- /dev/null +++ b/pvfactors/timeseries.py @@ -0,0 +1,487 @@ +# -*- coding: utf-8 -*- + +from pvlib import atmosphere, irradiance +from pvlib.tools import cosd, sind +from pvlib.irradiance import aoi_projection +import numpy as np +import pandas as pd +from pvfactors import (logging, print_progress) +from pvfactors.pvarray import Array +from multiprocessing import Pool, cpu_count +import time +from copy import deepcopy + + +LOGGER = logging.getLogger(__name__) +LOGGER.setLevel(logging.INFO) + + +COLS_TO_SAVE = ['q0', 'qinc', 'circumsolar_term', 'horizon_term', + 'direct_term', 'irradiance_term', 'isotropic_term', + 'reflection_term', 'horizon_band_shading_pct'] +DISTANCE_TOLERANCE = 1e-8 +idx_slice = pd.IndexSlice + + +def array_timeseries_calculate( + pvarray_parameters, timestamps, solar_zenith, solar_azimuth, array_tilt, + array_azimuth, dni, luminance_isotropic, luminance_circumsolar, + poa_horizon, poa_circumsolar): + """ + Calculate the view factor radiosity and irradiance terms for multiple + times. + The function inputs assume a diffuse sky dome as represented in the Perez + diffuse sky transposition model (from ``pvlib-python`` implementation). + + :param dict pvarray_parameters: parameters used to instantiate + ``pvarray.Array`` class + :param array-like timestamps: simulation timestamps + :param array-like solar_zenith: solar zenith angles + :param array-like solar_azimuth: solar azimuth angles + :param array-like array_tilt: pv module tilt angles + :param array-like array_azimuth: pv array azimuth angles + :param array-like dni: values for direct normal irradiance + :param array-like luminance_isotropic: luminance of the isotropic sky dome + :param array-like luminance_circumsolar: luminance of circumsolar area + :param array-like poa_horizon: POA irradiance horizon component + :param array-like poa_circumsolar: POA irradiance circumsolar component + + :return: ``df_registries``. + Concatenated form of the ``pvarray.Array`` surface registries. + :rtype: :class:`pandas.DataFrame` + """ + # Instantiate array + array = Array(**pvarray_parameters) + # We want to save the whole registry for each timestamp + list_registries = [] + # Use for printing progress + n = len(timestamps) + i = 1 + for idx, ts in enumerate(timestamps): + try: + if ((isinstance(solar_zenith[idx], float)) + & (solar_zenith[idx] <= 90.)): + # Run calculation only if daytime + array.calculate_radiosities_perez( + solar_zenith[idx], solar_azimuth[idx], array_tilt[idx], + array_azimuth[idx], dni[idx], luminance_isotropic[idx], + luminance_circumsolar[idx], poa_horizon[idx], + poa_circumsolar[idx]) + + # Save the whole registry + registry = deepcopy(array.surface_registry) + registry['timestamps'] = ts + list_registries.append(registry) + + except Exception as err: + LOGGER.debug("Unexpected error: {0}".format(err)) + + # Printing progress + print_progress(i, n, prefix='Progress:', suffix='Complete', + bar_length=50) + i += 1 + + # Concatenate all surface registries into one dataframe + if list_registries: + df_registries = pd.concat(list_registries, axis=0, join='outer') + else: + df_registries = None + + return df_registries + + +def perez_diffuse_luminance(timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, dni, dhi): + """ + Function used to calculate the luminance and the view factor terms from the + Perez diffuse light transposition model, as implemented in the + ``pvlib-python`` library. + This function was custom made to allow the calculation of the circumsolar + component on the back surface as well. Otherwise, the ``pvlib`` + implementation would ignore it. + + :param array-like timestamps: simulation timestamps + :param array-like array_tilt: pv module tilt angles + :param array-like array_azimuth: pv array azimuth angles + :param array-like solar_zenith: solar zenith angles + :param array-like solar_azimuth: solar azimuth angles + :param array-like dni: values for direct normal irradiance + :param array-like dhi: values for diffuse horizontal irradiance + :return: ``df_inputs``, dataframe with the following columns: + ['solar_zenith', 'solar_azimuth', 'array_tilt', 'array_azimuth', 'dhi', + 'dni', 'vf_horizon', 'vf_circumsolar', 'vf_isotropic', + 'luminance_horizon', 'luminance_circumsolar', 'luminance_isotropic', + 'poa_isotropic', 'poa_circumsolar', 'poa_horizon', 'poa_total_diffuse'] + :rtype: class:`pandas.DataFrame` + """ + # Create a dataframe to help filtering on all arrays + df_inputs = pd.DataFrame( + {'array_tilt': array_tilt, 'array_azimuth': array_azimuth, + 'solar_zenith': solar_zenith, 'solar_azimuth': solar_azimuth, + 'dni': dni, 'dhi': dhi}, + index=pd.DatetimeIndex(timestamps)) + + dni_et = irradiance.extraradiation(df_inputs.index.dayofyear) + am = atmosphere.relativeairmass(df_inputs.solar_zenith) + + # Need to treat the case when the sun is hitting the back surface of pvrow + aoi_proj = aoi_projection(df_inputs.array_tilt, df_inputs.array_azimuth, + df_inputs.solar_zenith, df_inputs.solar_azimuth) + sun_hitting_back_surface = ((aoi_proj < 0) & + (df_inputs.solar_zenith <= 90)) + df_inputs_back_surface = df_inputs.loc[sun_hitting_back_surface] + # Reverse the surface normal to switch to back-surface circumsolar calc + df_inputs_back_surface.loc[:, 'array_azimuth'] -= 180. + df_inputs_back_surface.loc[:, 'array_azimuth'] = np.mod( + df_inputs_back_surface.loc[:, 'array_azimuth'], 360. + ) + df_inputs_back_surface.loc[:, 'array_tilt'] = ( + 180. - df_inputs_back_surface.array_tilt) + + if df_inputs_back_surface.shape[0] > 0: + # Use recursion to calculate circumsolar luminance for back surface + df_inputs_back_surface = perez_diffuse_luminance( + *breakup_df_inputs(df_inputs_back_surface)) + + # Calculate Perez diffuse components + diffuse_poa, components = irradiance.perez(df_inputs.array_tilt, + df_inputs.array_azimuth, + df_inputs.dhi, df_inputs.dni, + dni_et, + df_inputs.solar_zenith, + df_inputs.solar_azimuth, + am, + return_components=True) + + # Calculate Perez view factors: + a = aoi_projection(df_inputs.array_tilt, df_inputs.array_azimuth, + df_inputs.solar_zenith, df_inputs.solar_azimuth) + a = np.maximum(a, 0) + b = cosd(df_inputs.solar_zenith) + b = np.maximum(b, cosd(85)) + + vf_perez = pd.DataFrame( + np.array([ + sind(df_inputs.array_tilt), + a / b, + (1. + cosd(df_inputs.array_tilt)) / 2. + ]).T, + index=df_inputs.index, + columns=['vf_horizon', 'vf_circumsolar', 'vf_isotropic'] + ) + + # Calculate diffuse luminance + luminance = pd.DataFrame( + np.array([ + components['horizon'] / vf_perez['vf_horizon'], + components['circumsolar'] / vf_perez['vf_circumsolar'], + components['isotropic'] / vf_perez['vf_isotropic'] + ]).T, + index=df_inputs.index, + columns=['luminance_horizon', 'luminance_circumsolar', + 'luminance_isotropic'] + ) + luminance.loc[diffuse_poa == 0, :] = 0. + + # Format components column names + components = components.rename(columns={'isotropic': 'poa_isotropic', + 'circumsolar': 'poa_circumsolar', + 'horizon': 'poa_horizon'}) + + df_inputs = pd.concat([df_inputs, components, vf_perez, luminance, + diffuse_poa], + axis=1, join='outer') + df_inputs = df_inputs.rename(columns={0: 'poa_total_diffuse'}) + + # Adjust the circumsolar luminance when it hits the back surface + if df_inputs_back_surface.shape[0] > 0: + df_inputs.loc[sun_hitting_back_surface, 'luminance_circumsolar'] = ( + df_inputs_back_surface.loc[:, 'luminance_circumsolar'] + ) + return df_inputs + + +def calculate_custom_perez_transposition(timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, + dni, dhi): + """ + Calculate custom perez transposition: some pre-processing is necessary in + order to get the circumsolar component when the sun is hitting the back + surface as well. + + :param array-like timestamps: simulation timestamps + :param array-like array_tilt: pv module tilt angles + :param array-like array_azimuth: pv array azimuth angles + :param array-like solar_zenith: solar zenith angles + :param array-like solar_azimuth: solar azimuth angles + :param array-like dni: values for direct normal irradiance + :param array-like dhi: values for diffuse horizontal irradiance + :return: ``df_custom_perez``, dataframe with the following columns: + ['solar_zenith', 'solar_azimuth', 'array_tilt', 'array_azimuth', 'dhi', + 'dni', 'vf_horizon', 'vf_circumsolar', 'vf_isotropic', + 'luminance_horizon', 'luminance_circumsolar', 'luminance_isotropic', + 'poa_isotropic', 'poa_circumsolar', 'poa_horizon', 'poa_total_diffuse'] + :rtype: class:`pandas.DataFrame` + """ + # Pre-process df_inputs to use the expected format of pvlib's perez model: + # only positive tilt angles, and switching azimuth angles + array_azimuth_processed = deepcopy(array_azimuth) + array_azimuth_processed[array_tilt < 0] = np.remainder( + array_azimuth[array_tilt < 0] + 180., + 360.) + array_tilt_processed = np.abs(array_tilt) + + # Calculate the perez inputs + df_custom_perez = perez_diffuse_luminance(timestamps, array_tilt_processed, + array_azimuth_processed, + solar_zenith, solar_azimuth, dni, dhi) + + return df_custom_perez + + +def calculate_radiosities_serially_perez(args): + """ Calculate timeseries results of simulation: run both custom Perez + diffuse light transposition calculations and ``pvarray.Array`` timeseries + calculation. + + :param args: tuple of arguments used to run the timeseries calculation. + List in order: ``pvarray_parameters``, ``timestamps``, + ``solar_zenith``, ``solar_azimuth``, ``array_tilt``, ``array_azimuth``, + ``dni``, ``dhi``. + All 1-dimensional arrays. + :return: ``df_registries``, ``df_custom_perez``; dataframes containing + the concatenated and timestamped ``pvarray.Array.surface_registry`` + values, and the custom Perez inputs used for it. + :rtype: both class:`pandas.DataFrame` + + """ + # Get arguments + (pvarray_parameters, timestamps, solar_zenith, solar_azimuth, + array_tilt, array_azimuth, dni, dhi) = args + + # Run custom perez transposition: in order to get circumsolar on back + # surface too + df_custom_perez = calculate_custom_perez_transposition( + timestamps, array_tilt, array_azimuth, solar_zenith, solar_azimuth, + dni, dhi) + + # Get the necessary inputs + luminance_isotropic = df_custom_perez.luminance_isotropic.values + luminance_circumsolar = df_custom_perez.luminance_circumsolar.values + poa_horizon = df_custom_perez.poa_horizon.values + poa_circumsolar = df_custom_perez.poa_circumsolar.values + + # Run timeseries calculation + df_registries = array_timeseries_calculate( + pvarray_parameters, timestamps, solar_zenith, solar_azimuth, + array_tilt, array_azimuth, dni, luminance_isotropic, + luminance_circumsolar, poa_horizon, poa_circumsolar) + + return df_registries, df_custom_perez + + +def calculate_radiosities_parallel_perez( + pvarray_parameters, timestamps, solar_zenith, solar_azimuth, + array_tilt, array_azimuth, dni, dhi, n_processes=None): + """ Calculate timeseries results of simulation in parallel: + run both custom Perez diffuse light transposition calculations and + ``pvarray.Array`` timeseries calculation. + + :param dict pvarray_parameters: parameters used to instantiate + ``pvarray.Array`` class + :param array-like timestamps: simulation timestamps + :param array-like solar_zenith: solar zenith angles + :param array-like solar_azimuth: solar azimuth angles + :param array-like array_tilt: pv module tilt angles + :param array-like array_azimuth: pv array azimuth angles + :param array-like dni: values for direct normal irradiance + :param array-like dhi: values for diffuse horizontal irradiance + :param int n_processes: (optional, default ``None`` = use all) number of + processes to use. Default value will be the total number of processors + on the machine. + :return: ``df_registries``, ``df_custom_perez``; dataframes containing + the concatenated and timestamped ``pvarray.Array.surface_registry`` + values, and the custom Perez inputs used for it. + :rtype: both class:`pandas.DataFrame` + """ + + # Choose number of workers + if n_processes is None: + n_processes = cpu_count() + + # Split all arguments according to number of processes + (list_timestamps, list_array_azimuth, list_array_tilt, + list_solar_zenith, list_solar_azimuth, list_dni, list_dhi) = map( + np.array_split, + [timestamps, array_azimuth, array_tilt, + solar_zenith, solar_azimuth, dni, dhi], + [n_processes] * 7) + list_parameters = [pvarray_parameters] * n_processes + # Zip all the arguments together + list_args = zip(*(list_parameters, list_timestamps, list_solar_zenith, + list_solar_azimuth, list_array_tilt, list_array_azimuth, + list_dni, list_dhi)) + + # Start multiprocessing + pool = Pool(n_processes) + start = time.time() + results = pool.map(calculate_radiosities_serially_perez, list_args) + end = time.time() + pool.close() + pool.join() + + LOGGER.info("Parallel calculation elapsed time: %s sec" % str(end - start)) + + results_grouped = zip(*results) + results_concat = map(lambda x: pd.concat(x, axis=0, join='outer') + if x[0] is not None else None, + results_grouped) + + return results_concat + + +def get_average_pvrow_outputs(df_registries, values=COLS_TO_SAVE, + include_shading=True): + """ For each pvrow surface (front and back), calculate the weighted + average irradiance and shading values (weighted by length of surface). + + :param df_registries: timestamped and concatenated :class:`pvarray.Array` + surface registries + :type df_registries: ``pandas.DataFrame`` + :param list values: list of column names from the surface registries to + keep (optional, default=``COLS_TO_SAVE``) + :param bool include_shading: flag to decide whether to keep column + indicating if surface has direct shading or not + :return: ``df_outputs``, dataframe with averaged values and flags + :rtype: ``pandas.DataFrame`` + """ + + weight_col = ['area'] + shade_col = ['shaded'] + indexes = ['timestamps', 'pvrow_index', 'surface_side'] + if include_shading: + final_cols = values + shade_col + else: + final_cols = values + + # Format registry to get averaged outputs for each surface type + df_outputs = ( + deepcopy(df_registries) + .query('line_type == "pvrow"') + .assign(pvrow_index=lambda x: x['pvrow_index'].astype(int)) + .loc[:, values + indexes + weight_col + shade_col] + # Weighted averages: make sure to close the col value in lambdas + .assign(**{col: lambda x, y=col: x[y] * x[weight_col[0]] + for col in values}) + .assign(shaded=lambda x: pd.to_numeric(x.shaded)) + .groupby(indexes) + .sum() + .assign(**{col: lambda x, y=col: x[y] / x[weight_col[0]] + for col in values}) + # summed up bool values, so there is shading if the shaded col > 0 + .assign(**{shade_col[0]: lambda x: x[shade_col[0]] > 0}) + # Now pivot data to the right format + .reset_index() + .melt(id_vars=indexes, value_vars=final_cols, var_name='term') + .pivot_table(index=['timestamps'], + columns=['pvrow_index', 'surface_side', 'term'], + values='value', + # The following works because there is no actual + # aggregation happening + aggfunc='first' # necessary to keep 'shaded' bool values + ) + ) + # Make sure numerical types are as expected + df_outputs.loc[:, idx_slice[:, :, values]] = df_outputs.loc[ + :, idx_slice[:, :, values]].astype(float) + + return df_outputs + + +def get_pvrow_segment_outputs(df_registries, values=COLS_TO_SAVE, + include_shading=True): + """For each discretized segment of a pvrow surface (front and back), + calculate the weighted average irradiance and shading values + (weighted by length of surface). + + :param df_registries: timestamped and concatenated :class:`pvarray.Array` + surface registries + :type df_registries: ``pandas.DataFrame`` + :param list values: list of column names from the surface registries to + keep (optional, default=``COLS_TO_SAVE``) + :param bool include_shading: flag to decide whether to keep column + indicating if surface has direct shading or not + :return: ``df_segments``, dataframe with averaged values and flags + :rtype: ``pandas.DataFrame`` + """ + + weight_col = ['area'] + shade_col = ['shaded'] + indexes = ['timestamps', 'pvrow_index', 'surface_side', + 'pvrow_segment_index'] + if include_shading: + final_cols = values + shade_col + else: + final_cols = values + + # Format registry to get averaged outputs for each surface type + df_segments = ( + deepcopy(df_registries) + .loc[~ np.isnan(df_registries.pvrow_segment_index).values, + values + indexes + weight_col + shade_col] + .assign( + pvrow_segment_index=lambda x: x['pvrow_segment_index'].astype(int), + pvrow_index=lambda x: x['pvrow_index'].astype(int), + shaded=lambda x: pd.to_numeric(x.shaded)) + # Weighted averages: make sure to close the col value in lambdas + # Include shaded and non shaded segments + .assign(**{col: lambda x, y=col: x[y] * x[weight_col[0]] + for col in values}) + .groupby(indexes) + .sum() + .assign(**{col: lambda x, y=col: x[y] / x[weight_col[0]] + for col in values}) + # summed up bool values, so there is shading if the shaded col > 0 + .assign(**{shade_col[0]: lambda x: x[shade_col[0]] > 0}) + # Now pivot data to the right format + .reset_index() + .melt(id_vars=indexes, value_vars=final_cols, var_name='term') + .pivot_table(index=['timestamps'], + columns=['pvrow_index', 'surface_side', + 'pvrow_segment_index', 'term'], + values='value', + # The following works because there is no actual + # aggregation happening + aggfunc='first' # necessary to keep 'shaded' bool values + ) + ) + # Make sure numerical types are as expected + df_segments.loc[:, idx_slice[:, :, :, values]] = df_segments.loc[ + :, idx_slice[:, :, :, values]].astype(float) + + return df_segments + + +def breakup_df_inputs(df_inputs): + """ It is sometimes easier to provide a dataframe than a list of arrays: + this function does the job of breaking up the dataframe into a list of + expected 1-dim arrays + + :param df_inputs: timestamp-indexed dataframe with following columns: + 'array_azimuth', 'array_tilt', 'solar_zenith', 'solar_azimuth', + 'dni', 'dhi' + :type df_inputs: ``pandas.DataFrame`` + :return: ``timestamps``, ``array_tilt``, ``array_azimuth``, + ``solar_zenith``, ``solar_azimuth``, ``dni``, ``dhi`` + :rtype: all 1-dim arrays + """ + timestamps = pd.to_datetime(df_inputs.index) + array_azimuth = df_inputs.array_azimuth.values + array_tilt = df_inputs.array_tilt.values + solar_zenith = df_inputs.solar_zenith.values + solar_azimuth = df_inputs.solar_azimuth.values + dni = df_inputs.dni.values + dhi = df_inputs.dhi.values + + return (timestamps, array_tilt, array_azimuth, + solar_zenith, solar_azimuth, dni, dhi) diff --git a/pvfactors/tools.py b/pvfactors/tools.py deleted file mode 100644 index 2abe6266..00000000 --- a/pvfactors/tools.py +++ /dev/null @@ -1,603 +0,0 @@ -# -*- coding: utf-8 -*- - -import sys -from pvlib import atmosphere, irradiance -from pvlib.tools import cosd, sind -from pvlib.irradiance import aoi_projection -import numpy as np -import pandas as pd -from pvfactors import (logging, PVFactorsError) -from pvfactors.pvarray import Array -from multiprocessing import Pool, cpu_count -import time -from copy import deepcopy - -DISTANCE_TOLERANCE = 1e-8 - -LOGGER = logging.getLogger(__name__) -LOGGER.setLevel(logging.INFO) - - -# Define colors used for plotting the 2D arrays -COLOR_dic = { - 'i': '#FFBB33', - 's': '#A7A49D', - 't': '#6699cc', - 'pvrow_illum': '#6699cc', - 'pvrow_shaded': '#ff0000', - 'ground_shaded': '#A7A49D', - 'ground_illum': '#FFBB33' -} - -COLS_TO_SAVE = ['q0', 'qinc', 'circumsolar_term', 'horizon_term', - 'direct_term', 'irradiance_term', 'isotropic_term', - 'reflection_term', 'horizon_band_shading_pct'] - -idx_slice = pd.IndexSlice - - -def plot_pvarray(ax, pvarray, line_types_selected=None, fontsize=20): - """ - Plot a :class:`pvarray.Array` object's shapely geometries based on its - :attr:`pvarray.Array.surface_registry`. The difference with - :func:`tools.plot_line_registry` is that here the user will see the - differences between PV row sides, the discretized elements, etc. - - :param ax: :class:`matplotlib.axes.Axes` object to use for the plot - :param pvarray: :class:`pvarray.Array` object to plot - :param list line_types_selected: parameter used to select a subset of - 'line_type' to plot; e.g. 'pvrow' or 'ground' - :return: None (``ax`` is updated) - """ - - # FIXME: repeating code from plot_line_registry - surface_registry = pvarray.surface_registry.copy() - plot_array_from_registry(ax, surface_registry, - line_types_selected=line_types_selected) - - # Plot details - distance = pvarray.pvrow_distance - height = pvarray.pvrow_height - n_pvrows = pvarray.n_pvrows - ax.set_xlim(- 0.5 * distance, (n_pvrows - 0.5) * distance) - ax.set_ylim(-height, 2 * height) - ax.set_title("PV Array", fontsize=fontsize) - - -def plot_array_from_registry(ax, registry, line_types_selected=None, - fontsize=20): - """ - Plot a 2D PV array geometry based using a registry input. - - :param matplotlib.axes.Axes ax: axes to use for the plot - :param pd.DataFrame registry: :class:`pvarray.Array` object to plot - :param list line_types_selected: parameter used to select a subset of - 'line_type' to plot; e.g. 'pvrow' or 'ground' - :return: None (``ax`` is updated) - """ - - registry = registry.copy() - registry.loc[:, 'color'] = ( - registry.line_type.values + '_' - + np.where(registry.shaded.values, 'shaded', 'illum')) - # TODO: distance may not exist - if line_types_selected: - for line_type in line_types_selected: - surface_reg_selected = registry.loc[ - registry.line_type == line_type, :] - for index, row in surface_reg_selected.iterrows(): - LOGGER.debug("Plotting %s", row['line_type']) - plot_coords(ax, row['geometry']) - plot_bounds(ax, row['geometry']) - plot_line(ax, row['geometry'], row['style'], - row['shading_type']) - else: - for index, row in registry.iterrows(): - LOGGER.debug("Plotting %s", row['line_type']) - plot_coords(ax, row['geometry']) - plot_bounds(ax, row['geometry']) - plot_line(ax, row['geometry'], row['style'], row['color']) - - ax.axis('equal') - ax.set_xlabel("x [m]", fontsize=fontsize) - ax.set_ylabel("y [m]", fontsize=fontsize) - - -def perez_diffuse_luminance(df_inputs): - """ - Function used to calculate the luminance and the view factor terms from the - Perez diffuse light transposition model, as implemented in the - ``pvlib-python`` library. - - :param df_inputs: class:`pandas.DataFrame` with following columns: - ['solar_zenith', 'solar_azimuth', 'array_tilt', 'array_azimuth', 'dhi', - 'dni']. Units are: ['deg', 'deg', 'deg', 'deg', 'W/m2', 'W/m2'] - :return: class:`pandas.DataFrame` with the following columns: - ['solar_zenith', 'solar_azimuth', 'array_tilt', 'array_azimuth', 'dhi', - 'dni', 'vf_horizon', 'vf_circumsolar', 'vf_isotropic', - 'luminance_horizon', 'luminance_circumsolar', 'luminance_isotropic', - 'poa_isotropic', 'poa_circumsolar', 'poa_horizon', 'poa_total_diffuse'] - """ - - dni_et = irradiance.extraradiation(df_inputs.index.dayofyear) - am = atmosphere.relativeairmass(df_inputs.solar_zenith) - - # Need to treat the case when the sun is hitting the back surface of pvrow - aoi_proj = aoi_projection(df_inputs.array_tilt, df_inputs.array_azimuth, - df_inputs.solar_zenith, df_inputs.solar_azimuth) - sun_hitting_back_surface = ((aoi_proj < 0) & - (df_inputs.solar_zenith <= 90)) - df_inputs_back_surface = df_inputs.loc[sun_hitting_back_surface] - # Reverse the surface normal to switch to back-surface circumsolar calc - df_inputs_back_surface.loc[:, 'array_azimuth'] -= 180. - df_inputs_back_surface.loc[:, 'array_azimuth'] = np.mod( - df_inputs_back_surface.loc[:, 'array_azimuth'], 360. - ) - df_inputs_back_surface.loc[:, 'array_tilt'] = ( - 180. - df_inputs_back_surface.array_tilt) - - if df_inputs_back_surface.shape[0] > 0: - # Use recursion to calculate circumsolar luminance for back surface - df_inputs_back_surface = perez_diffuse_luminance( - df_inputs_back_surface) - - # Calculate Perez diffuse components - diffuse_poa, components = irradiance.perez(df_inputs.array_tilt, - df_inputs.array_azimuth, - df_inputs.dhi, df_inputs.dni, - dni_et, - df_inputs.solar_zenith, - df_inputs.solar_azimuth, - am, - return_components=True) - - # Calculate Perez view factors: - a = aoi_projection(df_inputs.array_tilt, df_inputs.array_azimuth, - df_inputs.solar_zenith, df_inputs.solar_azimuth) - a = np.maximum(a, 0) - b = cosd(df_inputs.solar_zenith) - b = np.maximum(b, cosd(85)) - - vf_perez = pd.DataFrame( - np.array([ - sind(df_inputs.array_tilt), - a / b, - (1. + cosd(df_inputs.array_tilt)) / 2. - ]).T, - index=df_inputs.index, - columns=['vf_horizon', 'vf_circumsolar', 'vf_isotropic'] - ) - - # Calculate diffuse luminance - luminance = pd.DataFrame( - np.array([ - components['horizon'] / vf_perez['vf_horizon'], - components['circumsolar'] / vf_perez['vf_circumsolar'], - components['isotropic'] / vf_perez['vf_isotropic'] - ]).T, - index=df_inputs.index, - columns=['luminance_horizon', 'luminance_circumsolar', - 'luminance_isotropic'] - ) - luminance.loc[diffuse_poa == 0, :] = 0. - - # Format components column names - components = components.rename(columns={'isotropic': 'poa_isotropic', - 'circumsolar': 'poa_circumsolar', - 'horizon': 'poa_horizon'}) - - df_inputs = pd.concat([df_inputs, components, vf_perez, luminance, - diffuse_poa], - axis=1, join='outer') - df_inputs = df_inputs.rename(columns={0: 'poa_total_diffuse'}) - - # Adjust the circumsolar luminance when it hits the back surface - if df_inputs_back_surface.shape[0] > 0: - df_inputs.loc[sun_hitting_back_surface, 'luminance_circumsolar'] = ( - df_inputs_back_surface.loc[:, 'luminance_circumsolar'] - ) - return df_inputs - - -def calculate_radiosities_serially_simple(array, df_inputs): - """ - /!\ DEPRECATED - ============== - - Calculate the view factor radiosity and irradiance terms for multiple times. - The calculations will be sequential, and they will assume a completely - isotropic sky dome. - - :param array: :class:`pvarray.Array` object already configured and - instantiated - :param df_inputs: :class:`pandas.DataFrame` with following columns: - ['solar_zenith', 'solar_azimuth', 'array_tilt', 'array_azimuth', 'dhi', - 'dni']. Units are: ['deg', 'deg', 'deg', 'deg', 'W/m2', 'W/m2'] - :return: ``df_outputs, df_bifacial_gain``; :class:`pandas.DataFrame` objects - where ``df_outputs`` contains *averaged* irradiance terms for all PV row - sides and at each time stamp; ``df_bifacial_gain`` contains the - calculation of back-surface over front-surface irradiance for all PV - rows and at each time stamp. - """ - LOGGER.warning("``calculate_radiosities_serially_simple`` is deprecated") - # Create index df_outputs - iterables = [ - range(array.n_pvrows), - ['front', 'back'], - ['q0', 'qinc'] - ] - multiindex = pd.MultiIndex.from_product(iterables, - names=['pvrow', 'side', 'term']) - - df_outputs = pd.DataFrame(np.nan, columns=df_inputs.index, - index=multiindex) - df_outputs.sort_index(inplace=True) - df_outputs.loc['array_is_shaded', :] = np.nan - idx_slice = pd.IndexSlice - - n = df_inputs.shape[0] - i = 1 - - for idx, row in df_inputs.iterrows(): - try: - - if ((isinstance(row['solar_zenith'], float)) - & (row['solar_zenith'] <= 90.)): - array.calculate_radiosities_simple(row['solar_zenith'], - row['solar_azimuth'], - row['array_tilt'], - row['array_azimuth'], - row['dni'], row['dhi']) - - array.surface_registry.loc[:, 'q0'] = ( - array.surface_registry.loc[:, 'area'] - * array.surface_registry.q0) - array.surface_registry.loc[:, 'qinc'] = ( - array.surface_registry.loc[:, 'area'] - * array.surface_registry.qinc - ) - df_summed = array.surface_registry.groupby( - ['pvrow_index', 'surface_side']).sum() - df_avg_irradiance = df_summed.div(df_summed['area'], - axis=0).loc[ - idx_slice[:, :], ['q0', 'qinc']].sort_index() - df_outputs.loc[idx_slice[:, :, 'q0'], idx] = ( - df_avg_irradiance.loc[ - idx_slice[:, :], 'q0'].values - ) - df_outputs.loc[idx_slice[:, :, 'qinc'], - idx] = df_avg_irradiance.loc[ - idx_slice[:, :], 'qinc'].values - df_outputs.loc['array_is_shaded', idx] = ( - array.has_direct_shading) - except Exception as err: - LOGGER.debug("Unexpected error: {0}".format(err)) - - print_progress(i, n, prefix='Progress:', suffix='Complete', - bar_length=50) - i += 1 - - try: - bifacial_gains = (df_outputs.loc[ - idx_slice[:, 'back', 'qinc'], :].values - / df_outputs.loc[ - idx_slice[:, 'front', 'qinc'], :].values) - df_bifacial_gain = pd.DataFrame(bifacial_gains.T, - index=df_outputs.index, - columns=range(array.n_pvrows)) - except Exception as err: - LOGGER.warning("Error in calculation of bifacial gain %s" % err) - df_bifacial_gain = pd.DataFrame( - np.nan * np.ones((len(df_inputs.index), array.n_pvrows)), - index=df_inputs.index, - columns=range(array.n_pvrows)) - - return df_outputs, df_bifacial_gain - - -def calculate_radiosities_serially_perez(args): - """ - Calculate the view factor radiosity and irradiance terms for multiple times. - The calculations will be sequential, and they will assume a diffuse sky - dome as calculated in the Perez diffuse sky transposition model (from - ``pvlib-python`` implementation). - - :param args: tuple of at least two arguments: ``(arguments, df_inputs)``, - where ``arguments`` is a ``dict`` that contains the array parameters - used to instantiate a :class:`pvarray.Array` object, and ``df_inputs`` - is a :class:`pandas.DataFrame` with the following columns: - ['solar_zenith', 'solar_azimuth', 'array_tilt', 'array_azimuth', 'dhi', - 'dni'], and with the following units: ['deg', 'deg', 'deg', 'deg', - 'W/m2', 'W/m2']. A possible 3rd argument for the tuple is - ``save_segments``, which is a ``tuple`` of two elements used to save - all the irradiance terms calculated for one side of a PV row; e.g. - ``(1, 'front')`` the first element is an ``int`` for the PV row index, - and the second element a ``str`` to specify the side of the PV row, - 'front' or 'back' - :return: ``df_outputs, df_bifacial_gain, df_inputs_perez, ipoa_dict``; - :class:`pandas.DataFrame` objects and ``dict`` where ``df_outputs`` - contains *averaged* irradiance terms for all PV row sides and at each - time stamp; ``df_bifacial_gain`` contains the calculation of - back-surface over front-surface irradiance for all PV rows and at each - time stamp; ``df_inputs_perez`` contains the intermediate input and - output values from the Perez model calculation in ``pvlib-python``; - ``ipoa_dict`` is not ``None`` only when the ``save_segments`` input is - specified, and it is otherwise a ``dict`` where the keys are all the - calculated irradiance terms' names, and the values are - :class:`pandas.DataFrame` objects containing the calculated values for - all the segments of the PV string side (it is a way of getting detailed - values instead of averages) - """ - - if len(args) == 3: - arguments, df_inputs, save_segments = args - else: - arguments, df_inputs = args - save_segments = None - - array = Array(**arguments) - # Pre-process df_inputs to use the expected format of pvlib's perez model: - # only positive tilt angles, and switching azimuth angles - df_inputs_before_perez = df_inputs.copy() - df_inputs_before_perez.loc[df_inputs.array_tilt < 0, 'array_azimuth'] = ( - np.remainder( - df_inputs_before_perez.loc[df_inputs.array_tilt < 0, - 'array_azimuth'] + 180., - 360.) - ) - df_inputs_before_perez.array_tilt = np.abs(df_inputs_before_perez - .array_tilt) - - # Calculate the perez inputs - df_inputs_perez = perez_diffuse_luminance(df_inputs_before_perez) - - # Post process: in vf model tilt angles can be negative and azimuth is - # generally fixed - df_inputs_perez.loc[:, ['array_azimuth', 'array_tilt']] = ( - df_inputs.loc[:, ['array_azimuth', 'array_tilt']] - ) - - # We want to save the whole registry for each timestamp - list_registries = [] - - # Use for printing progress - n = df_inputs_perez.shape[0] - i = 1 - - for idx, row in df_inputs_perez.iterrows(): - try: - if ((isinstance(row['solar_zenith'], float)) - & (row['solar_zenith'] <= 90.)): - array.calculate_radiosities_perez(row['solar_zenith'], - row['solar_azimuth'], - row['array_tilt'], - row['array_azimuth'], - row['dni'], - row['luminance_isotropic'], - row['luminance_circumsolar'], - row['poa_horizon'], - row['poa_circumsolar']) - - # Save the whole registry - registry = deepcopy(array.surface_registry) - registry['timestamps'] = idx - list_registries.append(registry) - - except Exception as err: - LOGGER.debug("Unexpected error: {0}".format(err)) - - # Printing progress - print_progress(i, n, prefix='Progress:', suffix='Complete', - bar_length=50) - i += 1 - - # Save all the registries into 1 output - if list_registries: - df_registries = pd.concat(list_registries, axis=0, join='outer') - else: - df_registries = None - - return df_registries, df_inputs_perez - - -def calculate_radiosities_parallel_perez(array_parameters, df_inputs, - n_processes=None, save_segments=None): - """ - Calculate the view factor radiosity and irradiance terms for multiple times. - The calculations will be run in parallel on the different processors, and - they will assume a diffuse sky dome as calculated in the Perez diffuse - sky transposition model (from ``pvlib-python`` implementation). - This function uses :func:`tools.calculate_radiosities_serially_perez` - - :param dict array_parameters: contains the array parameters used to - instantiate a :class:`pvarray.Array` object - :param df_inputs: :class:`pandas.DataFrame` with the following columns: - ['solar_zenith', 'solar_azimuth', 'array_tilt', 'array_azimuth', 'dhi', - 'dni'], and with the following units: ['deg', 'deg', 'deg', 'deg', - 'W/m2', 'W/m2'] - :param int n_processes: number of processes to use. Default value will be - the total number of processors on the machine - :param tuple save_segments: ``tuple`` of two elements used to save all the - irradiance terms calculated for one side of a PV row; e.g. - ``(1, 'front')`` the first element is an ``int`` for the PV row index, - and the second element a ``str`` to specify the side of the PV row, - 'front' or 'back' - :return: concatenated outputs of - :func:`tools.calculate_radiosities_serially_perez` run and outputted in - the parallel processes - """ - - # Choose number of workers - if n_processes is None: - n_processes = cpu_count() - - # Define list of arguments for target function - list_df_inputs = np.array_split(df_inputs, n_processes, axis=0) - list_parameters = [array_parameters] * n_processes - list_save_segments = [save_segments] * n_processes - list_args = zip(*(list_parameters, list_df_inputs, list_save_segments)) - - # Start multiprocessing - pool = Pool(n_processes) - start = time.time() - results = pool.map(calculate_radiosities_serially_perez, list_args) - end = time.time() - pool.close() - pool.join() - - LOGGER.info("Parallel calculation elapsed time: %s sec" % str(end - start)) - - results_grouped = zip(*results) - results_concat = map(lambda x: pd.concat(x, axis=0, join='outer') - if x[0] is not None else None, - results_grouped) - - return results_concat - - -# Base functions used to plot the 2D array -def plot_coords(ax, ob): - try: - x, y = ob.xy - ax.plot(x, y, 'o', color='#999999', zorder=1) - except NotImplementedError: - for line in ob: - x, y = line.xy - ax.plot(x, y, 'o', color='#999999', zorder=1) - - -def plot_bounds(ax, ob): - # Check if shadow reduces to one point (for very specific sun alignment) - if len(ob.boundary) == 0: - x, y = ob.coords[0] - else: - x, y = zip(*list((p.x, p.y) for p in ob.boundary)) - ax.plot(x, y, 'o', color='#000000', zorder=1) - - -def plot_line(ax, ob, line_style, line_color): - try: - x, y = ob.xy - ax.plot(x, y, color=COLOR_dic[line_color], ls=line_style, alpha=0.7, - linewidth=3, solid_capstyle='round', zorder=2) - except NotImplementedError: - for line in ob: - x, y = line.xy - ax.plot(x, y, color=COLOR_dic[line_color], ls=line_style, alpha=0.7, - linewidth=3, solid_capstyle='round', zorder=2) - - -# Define function used for progress bar when running long simulations -# Borrowed from: https://gist.github.com/aubricus -def print_progress(iteration, total, prefix='', suffix='', decimals=1, - bar_length=100): - """ - Call in a loop to create terminal progress bar - @params: - iteration - Required : current iteration (Int) - total - Required : total iterations (Int) - prefix - Optional : prefix string (Str) - suffix - Optional : suffix string (Str) - decimals - Optional : positive number of decimals in percent - complete (Int) - bar_length - Optional : character length of bar (Int) - """ - format_str = "{0:." + str(decimals) + "f}" - percents = format_str.format(100 * (iteration / float(total))) - filled_length = int(round(bar_length * iteration / float(total))) - bar = '█' * filled_length + '-' * (bar_length - filled_length) - sys.stdout.write('\r%s |%s| %s%s %s' % (prefix, bar, percents, '%', - suffix)), - sys.stdout.flush() - if iteration == total: - sys.stdout.write('\n') - sys.stdout.flush() - - -def get_average_pvrow_outputs(df_registries, values=COLS_TO_SAVE): - """ Calculate surface side irradiance averages for the pvrows """ - - weight_col = ['area'] - shade_col = ['shaded'] - indexes = ['timestamps', 'pvrow_index', 'surface_side'] - - # Format registry to get averaged outputs for each surface type - df_outputs = ( - deepcopy(df_registries) - .query('line_type == "pvrow"') - .assign(pvrow_index=lambda x: x['pvrow_index'].astype(int)) - .loc[:, values + indexes + weight_col + shade_col] - # Calculate weighted averages: make sure to close the col value in lambdas - .assign(**{col: lambda x, y=col: x[y] * x[weight_col[0]] - for col in values}) - .assign(shaded=lambda x: pd.to_numeric(x.shaded)) - .groupby(indexes) - .sum() - .assign(**{col: lambda x, y=col: x[y] / x[weight_col[0]] - for col in values}) - # summed up bool values, so there is shading if the shaded col > 0 - .assign(**{shade_col[0]: lambda x: x[shade_col[0]] > 0}) - # Now pivot data to the right format - .reset_index() - .melt(id_vars=indexes + shade_col, value_vars=values, var_name='term') - .pivot_table(index=['timestamps'] + shade_col, - columns=['pvrow_index', 'surface_side', 'term'], - values='value') - # Get shading as a column - .reset_index(level=1) - ) - - return df_outputs - - -def get_bifacial_gain_outputs(df_outputs): - """ Calculate irradiance bifacial gain for all pvrows """ - pass - - -def get_pvrow_segment_outputs(df_registries, values=COLS_TO_SAVE, - include_shading=True): - """ Get only pvrow segment outputs """ - - weight_col = ['area'] - shade_col = ['shaded'] - indexes = ['timestamps', 'pvrow_index', 'surface_side', 'pvrow_segment_index'] - if include_shading: - final_cols = values + shade_col - else: - final_cols = values - - # Format registry to get averaged outputs for each surface type - df_segments = ( - deepcopy(df_registries) - .loc[~ np.isnan(df_registries.pvrow_segment_index).values, - values + indexes + weight_col + shade_col] - .assign(pvrow_segment_index=lambda x: x['pvrow_segment_index'].astype(int), - pvrow_index=lambda x: x['pvrow_index'].astype(int), - shaded=lambda x: pd.to_numeric(x.shaded)) - # Calculate weighted averages: make sure to close the col value in lambdas - # Include shaded and non shaded segments - .assign(**{col: lambda x, y=col: x[y] * x[weight_col[0]] - for col in values}) - .groupby(indexes) - .sum() - .assign(**{col: lambda x, y=col: x[y] / x[weight_col[0]] - for col in values}) - # summed up bool values, so there is shading if the shaded col > 0 - .assign(**{shade_col[0]: lambda x: x[shade_col[0]] > 0}) - # Now pivot data to the right format - .reset_index() - .melt(id_vars=indexes, value_vars=final_cols, var_name='term') - .pivot_table(index=['timestamps'], - columns=['pvrow_index', 'surface_side', 'pvrow_segment_index', - 'term'], - values='value', - # The following works because there is no actual aggregation happening - aggfunc='first' # necessary to keep 'shaded' bool values - ) - ) - # Make sure numerical types are as expected - df_segments.loc[:, idx_slice[:, :, :, values]] = df_segments.loc[:, idx_slice[:, :, :, values] - ].astype(float) - - return df_segments diff --git a/setup.py b/setup.py index 4eee30c6..b0b73191 100644 --- a/setup.py +++ b/setup.py @@ -31,7 +31,7 @@ 'future>=0.16.0', 'six>=1.11.0'] -TESTS_REQUIRE = ['pytest>=3.2.1'] +TESTS_REQUIRE = ['pytest>=3.2.1', 'pytest-mock>=1.10.0'] setup(name=DISTNAME, description=DESCRIPTION,