Merge pull request #147 from nismod/p3

P3

energy_demand/basic/conversions.py

@@ -80,6 +80,23 @@ def kwh_to_gwh(kwh):
 
     return gwh
 
+def gwh_to_kwh(gwh):
+    """"Conversion of MW to GWh
+
+    Arguments
+    ---------
+    kwh : float
+        Kilowatthours
+
+    Return
+    ------
+    gwh : float
+        Gigawatthours
+    """
+    kwh = gwh *  1000000
+
+    return kwh
+
 def mw_to_gwh(megawatt, number_of_hours):
     """"Conversion of MW to GWh
 

energy_demand/main.py

@@ -1,7 +1,7 @@
 """Allows to run HIRE locally outside the SMIF framework
 
 write_national_results_amman
-
+move data config folder into data
 Note
 ----
 Always execute from root folder. (e.g. energy_demand/energy_demand/main.py
@@ -105,7 +105,8 @@ def energy_demand_model(
     user_defined_simulation_end_yr = config['CONFIG']['user_defined_simulation_end_yr']
 
     # Simulated yrs
-    sim_yrs = [base_yr, user_defined_simulation_end_yr]
+    #sim_yrs = [base_yr, user_defined_simulation_end_yr]
+    sim_yrs = [2015, 2025, 2035, 2045, 2050]
     weather_yr_scenario = 2015   # Default weather year
 
     if len(sys.argv) > 1: #user defined arguments are provide
@@ -292,7 +293,10 @@ def energy_demand_model(
         weather_realisation=weather_realisation,
         path_weather_data=path_weather_data,
         same_base_year_weather=False,
-        crit_temp_min_max=config['CRITERIA']['crit_temp_min_max'])
+        crit_temp_min_max=config['CRITERIA']['crit_temp_min_max'],
+        load_np=False,
+        load_parquet=False,
+        load_csv=True)
 
     #print(station_id_253)
     # Plot map with weather station
@@ -543,7 +547,7 @@ def energy_demand_model(
                     'fueltype_reg_peak_day')
 
                 # PLot only residential total regional annual demand and
-                write_data.write_residential_tot_demands(
+                '''write_data.write_residential_tot_demands(
                     sim_yr,
                     path_runs,
                     sim_obj.ed_residential_tot_reg_y,
@@ -553,7 +557,7 @@ def energy_demand_model(
                     "ed_fueltype_regs_yh",
                     path_runs,
                     sim_obj.ed_fueltype_regs_yh,
-                    "result_tot_submodels_fueltypes")
+                    "result_tot_submodels_fueltypes")'''
             else:
                 write_data.write_residential_tot_demands(
                     sim_yr,

energy_demand/plotting/fig_3_plot_over_time.py

@@ -0,0 +1,111 @@
+
+"""
+"""
+
+import os
+import logging
+import copy
+import math
+import numpy as np
+import geopandas as gpd
+import pandas as pd
+import palettable
+import matplotlib.pyplot as plt
+import matplotlib.patches as mpatches
+from matplotlib.patches import Circle
+from matplotlib.colors import LinearSegmentedColormap
+from collections import defaultdict
+
+from energy_demand.plotting import basic_plot_functions
+from energy_demand.plotting import fig_p2_weather_val, result_mapping
+from energy_demand.basic import basic_functions
+from energy_demand.technologies import tech_related
+from energy_demand.read_write import write_data
+from energy_demand.basic import conversions
+
+def scenario_over_time(
+        scenario_result_container,
+        sim_yrs,
+        fig_name,
+        result_path
+    ):
+
+    color_list = [
+        'red', 'green', 'orange', '#37AB65',
+        '#C0E4FF', '#3DF735', '#AD6D70', '#EC2504',
+        '#8C0B90', '#27B502', '#7C60A8', '#CF95D7', '#F6CC1D']
+
+    for cnt, i in enumerate(scenario_result_container):
+        scenario_name = i['scenario_name']
+        ed_reg_tot_y = i['national_peak']
+        # dataframe with national peak (columns= simulation year, row: Realisation) 
+
+        # Calculate quantiles
+        quantile_95 = 0.95
+        quantile_05 = 0.05
+
+        color = color_list[cnt]
+
+        # Calculate average across all weather scenarios
+        mean_ed_reg_tot_y = ed_reg_tot_y.mean(axis=0)
+
+        # Standard deviation over all realisations
+        df_q_05 = ed_reg_tot_y.quantile(quantile_05)
+        df_q_95 = ed_reg_tot_y.quantile(quantile_95)
+
+        # --------------------
+        # Try to smooth lines
+        # --------------------
+        try:
+            sim_yrs_smoothed, mean_ed_reg_tot_y_smoothed = basic_plot_functions.smooth_data(sim_yrs, mean_ed_reg_tot_y, num=40000)
+            _, df_q_05_smoothed = basic_plot_functions.smooth_data(sim_yrs, df_q_05, num=40000)
+            _, df_q_95_smoothed = basic_plot_functions.smooth_data(sim_yrs, df_q_95, num=40000)
+
+            mean_ed_reg_tot_y = pd.Series(mean_ed_reg_tot_y_smoothed, sim_yrs_smoothed)
+            sim_yrs = pd.Series(sim_yrs_smoothed, sim_yrs_smoothed)
+            #sim_yrs = list(sim_yrs_smoothed)
+            df_q_05 = pd.Series(df_q_05_smoothed, sim_yrs_smoothed)
+            df_q_95 = pd.Series(df_q_95_smoothed, sim_yrs_smoothed)
+        except:
+            pass
+
+
+        fig = plt.figure(figsize=basic_plot_functions.cm2inch(9, 8)) #width, height
+        ax = fig.add_subplot(1, 1, 1)
+
+        plt.plot(mean_ed_reg_tot_y, label=scenario_name, color=color)
+
+        # Plottin qunatilse and average scenario
+        df_q_05.plot.line(color=color, linestyle='--', linewidth=0.5, label="0.05")
+        df_q_95.plot.line(color=color, linestyle='--', linewidth=0.5, label="0.05")
+        plt.fill_between(
+            sim_yrs,
+            list(df_q_95),  #y1
+            list(df_q_05),  #y2
+            alpha=0.15,
+            facecolor=color,
+            label="uncertainty band")
+
+    #plt.ylim(0, y_lim_val)
+    plt.xlim(2015, 2050)
+
+    # --------
+    # Legend
+    # --------
+    legend = plt.legend(
+        title="tt",
+        prop={'size': 8},
+        loc='upper center',
+        bbox_to_anchor=(0.5, -0.05),
+        frameon=False)
+
+    # --------
+    # Labeling
+    # --------
+    plt.ylabel("peak")
+    plt.xlabel("year")
+    plt.title("tttt")
+
+    plt.tight_layout()
+    #plt.show()
+    plt.savefig(os.path.join(result_path, fig_name))

energy_demand/plotting/fig_3_weather_map.py

@@ -17,21 +17,25 @@
 from energy_demand.basic import basic_functions
 from energy_demand.technologies import tech_related
 from energy_demand.read_write import write_data
+from energy_demand.basic import conversions
 
 def total_annual_demand(
         df_data_input,
         path_shapefile_input,
         regions,
-        fueltypes_nr,
-        fueltypes,
         pop_data,
         simulation_yr_to_plot,
         result_path,
         fig_name,
+        field_to_plot,
+        unit='GW'
     ):
     """
     """
-    conversion_factor = 1000000 #  1000000 #GW to KW
+    if unit == 'GW':
+        conversion_factor = 1
+    if unit == 'kWh':
+        conversion_factor = conversions.gwh_to_kwh(gwh=1) #GW to KW
 
     df_data_input = df_data_input * conversion_factor
 
@@ -65,10 +69,12 @@ def total_annual_demand(
     print("TOT PERSON: " + str(tot_demand))
     print('AVERAGE KW per Person " '+ str(tot_demand / tot_person))
 
+    print(df_data_input)
     regional_statistics_columns = [
         'name',
         'mean',
-        'mean_norm_pop'] #,
+        'mean_norm_pop',
+        'std_dev']#
         #'diff_av_max',
         #'mean_pp',
         #'diff_av_max_pp',
@@ -82,7 +88,9 @@ def total_annual_demand(
         line_entry = [[
             str(region_name),
             mean[region_name],
-            mean_norm_pop[region_name], #diff_av_max,
+            mean_norm_pop[region_name],
+            std_dev[region_name]
+            #diff_av_max,
             #mean_peak_h_pp,
             #diff_av_max_pp,
             #std_dev_average_every_h,
@@ -113,11 +121,6 @@ def total_annual_demand(
     # Assign bin colors according to defined cmap and whether
     # plot with min_max values or only min/max values
     #bin_values = [0, 0.0025, 0.005, 0.0075, 0.01]
-
-    # Field to plot
-    field_to_plot = 'mean'
-    field_to_plot = 'mean_norm_pop'
-
     nr_of_intervals = 6
 
     bin_values = result_mapping.get_reasonable_bin_values_II(
@@ -148,23 +151,34 @@ def total_annual_demand(
         min_value,
         max_value)
 
-    plt.legend(
+    legend = plt.legend(
         handles=legend_handles,
-        title="tittel_elgend",
+        title="Unit: {} field: {}".format(unit, field_to_plot),
         prop={'size': 8},
         loc='upper center',
         bbox_to_anchor=(0.5, -0.05),
         frameon=False)
 
+    # Remove coordinates from figure
+    ax.set_yticklabels([])
+    ax.set_xticklabels([])
+
+    legend.get_title().set_fontsize(8)
+
     # PLot bins on plot
     '''plt.text(
         0,
         -20,
         bin_values[:-1],
         fontsize=8)'''
 
+    # --------
+    # Labeling
+    # --------
+    #plt.title("tttt")
+
     plt.tight_layout()
-    plt.show()
+    #plt.show()
 
-    plt.savefig(fig_name)
+    plt.savefig(os.path.join(result_path, fig_name))
     plt.close()

energy_demand/read_write/data_loader.py

@@ -901,7 +901,10 @@ def load_temp_data(
         weather_realisation,
         path_weather_data,
         same_base_year_weather=False,
-        crit_temp_min_max=False
+        crit_temp_min_max=False,
+        load_np=False,
+        load_parquet=False,
+        load_csv=True
     ):
     """Read in cleaned temperature and weather station data
 
@@ -933,10 +936,6 @@ def load_temp_data(
     """
     print("... loading temperatures", flush=True)
 
-    load_np = False
-    load_parquet = True
-    load_csv = False
-
     temp_data_short = defaultdict(dict)
     weather_stations_with_data = defaultdict(dict)
 

energy_demand/read_write/read_weather_results.py

@@ -9,6 +9,7 @@
 from energy_demand.basic import basic_functions
 from energy_demand.basic import lookup_tables
 from energy_demand.read_write import read_data
+from energy_demand.technologies import tech_related
 
 def read_in_weather_results(
         path_result,
@@ -29,7 +30,8 @@ def read_in_weather_results(
     """
     logging.info("... Reading in results")
 
-    lookups = lookup_tables.basic_lookups()
+    fueltype_int = tech_related.get_fueltype_int('electricity')
+
 
     results_container = {}
 
@@ -49,6 +51,20 @@ def read_in_weather_results(
         os.path.join('simulation_results', path_result), 'only_peak')
 
     print(results_container['ed_reg_peakday'][2015].shape)
+    results_container['ed_reg_peakday_peak_hour'] = {}
+    results_container['national_peak'] = {}
+    for year in results_container['ed_reg_peakday']:
+
+        # Get peak demand of each region
+        results_container['ed_reg_peakday_peak_hour'][year] = results_container['ed_reg_peakday'][year].max(axis=2)
+
+        # Get national peak
+        national_demand_per_hour = results_container['ed_reg_peakday'][year].sum(axis=1) #Aggregate houraly across all regions
+
+        # Get maximum hour for electricity demand
+        max_hour = national_demand_per_hour[fueltype_int].argmax()
+
+        results_container['national_peak'][year] = national_demand_per_hour[:, max_hour]
 
     logging.info("... Reading in results finished")
     return results_container