New_heatdem

.gitignore

@@ -2,7 +2,8 @@ RawDataProcessing/Data/*
 RawDataProcessing/Output/*
 !RawDataProcessing/Data/IncFilePreSuf/README.md
 !RawDataProcessing/Output/Figures/.keepoutfolder
-conversion_dictionaries.pkl
+exo_elec_dem_conversion_dictionaries.pkl
+exo_heat_dem_conversion_dictionaries.pkl
 .snakemake
 __pycache__
 .new_env/*

RawDataProcessing/Analysis/dot_source.txt

@@ -0,0 +1,24 @@
+digraph snakemake_dag {
+    graph[bgcolor=white, margin=0];
+    node[shape=box, style=rounded, fontname=sans,                 fontsize=10, penwidth=2];
+    edge[penwidth=2, color=grey];
+	0[label = "all", color = "0.33 0.6 0.85", style="rounded"];
+	1[label = "exo_electricity_demand", color = "0.40 0.6 0.85", style="rounded"];
+	2[label = "create_conversion_dictionaries", color = "0.53 0.6 0.85", style="rounded"];
+	3[label = "format_energinet_electricity_data", color = "0.00 0.6 0.85", style="rounded,dashed"];
+	4[label = "get_grid_from_Balmorel", color = "0.07 0.6 0.85", style="rounded"];
+	5[label = "format_dkstat_transport_data", color = "0.47 0.6 0.85", style="rounded,dashed"];
+	6[label = "exo_heat_demand", color = "0.13 0.6 0.85", style="rounded"];
+	7[label = "format_vpdk21_data", color = "0.60 0.6 0.85", style="rounded,dashed"];
+	8[label = "format_dkstat_industry_data", color = "0.27 0.6 0.85", style="rounded,dashed"];
+	1 -> 0
+	4 -> 0
+	5 -> 0
+	6 -> 0
+	2 -> 1
+	3 -> 1
+	2 -> 6
+	7 -> 6
+	8 -> 6
+	3 -> 6
+}            

RawDataProcessing/Analysis/plot_dag.py

@@ -5,23 +5,9 @@
 
 from graphviz import Source
 
-# Copy paste the output of snakemake --dag to within the dot_source below
-dot_source = '''
-digraph snakemake_dag {
-    graph[bgcolor=white, margin=0];
-    node[shape=box, style=rounded, fontname=sans,                 fontsize=10, penwidth=2];
-    edge[penwidth=2, color=grey];
-        0[label = "all", color = "0.56 0.6 0.85", style="rounded"];
-        1[label = "exo_electricity_demand", color = "0.44 0.6 0.85", style="rounded"];
-        2[label = "create_conversion_dictionaries", color = "0.11 0.6 0.85", style="rounded,dashed"];
-        3[label = "format_energinet_electricity_data_to_xarray", color = "0.33 0.6 0.85", style="rounded"];
-        4[label = "get_grid_from_Balmorel", color = "0.00 0.6 0.85", style="rounded"];
-        1 -> 0
-        4 -> 0
-        2 -> 1
-        3 -> 1
-}
-'''
+# Do snakemake --dag > analysis/dot_source.txt to within the dot_source below
+with open('Analysis/dot_source.txt', 'r') as f:
+    dot_source = f.read()
 
 # Create a Source object from the DOT source
 graph = Source(dot_source)

RawDataProcessing/Modules/Submodules/create_name_conversion_dict.py

@@ -3,8 +3,8 @@
 import pickle
 
 
-
-conversion_dictionaries = {
+# 1. Exogenous Electricity Demand xarray
+exo_elec_dem_conversion_dictionaries = {
 
 
     'coord_names' : {'municipality' : 'R',
@@ -33,5 +33,39 @@
 
 }
 
-with open('Modules/Submodules/conversion_dictionaries.pkl', 'wb') as f:
-    pickle.dump(conversion_dictionaries, f)
+with open('Modules/Submodules/exo_elec_dem_conversion_dictionaries.pkl', 'wb') as f:
+    pickle.dump(exo_elec_dem_conversion_dictionaries, f)
+
+
+# 2. Exogenous Heat Demand xarray
+exo_heat_dem_conversion_dictionaries = {
+
+
+    'coord_names' : {'municipality' : 'A',
+                                'user' : 'DHUSER',
+                                'year' : 'Y'},
+
+    'coord_element_names' : {'user' : {'district_heating' : 'RESH',
+                                       'individual' : 'RESIDENTIAL',
+                            'industry_phl' : 'IND-PHL',
+                            'industry_phm' : 'IND-PHM',
+                            'industry_phh' : 'IND-PHH'},
+                             'municipality' : {'æ' : 'ae',
+                                               'ø' : 'oe',
+                                               'å' : 'aa',
+                                               'Æ' : 'Ae',
+                                               'Ø' : 'Oe',
+                                               'Å' : 'Aa'
+                                 }},
+
+    'week_to_seasons' : ['S0%d'%i for i in range(1, 10)] +\
+                        ['S%d'%i for i in range(10, 53)],
+
+    'hour_to_terms'  : ['T00%d'%i for i in range(1, 10)] +\
+                       ['T0%d'%i for i in range(10, 100)] +\
+                       ['T%d'%i for i in range(100, 169)]
+
+}
+
+with open('Modules/Submodules/exo_heat_dem_conversion_dictionaries.pkl', 'wb') as f:
+    pickle.dump(exo_heat_dem_conversion_dictionaries, f)

RawDataProcessing/Modules/Submodules/utils.py

@@ -10,15 +10,16 @@
 
 import pandas as pd
 import xarray as xr
-from typing import Union
-
+from typing import Tuple
+import pickle
+import numpy as np
 
 #%% ------------------------------- ###
 ###     1. Conversion Functions     ###
 ### ------------------------------- ###
 
 # 1. Convert coord names and elements
-def convert_coordname_elements(dataset: Union[xr.Dataset, pd.DataFrame],
+def convert_coordname_elements(dataset: Tuple[xr.Dataset, pd.DataFrame],
                           dimension_map: dict,
                           element_map: dict,
                           print_before_and_after: bool = False
@@ -52,3 +53,51 @@ def convert_coordname_elements(dataset: Union[xr.Dataset, pd.DataFrame],
 
     return new_dataset
 
+
+def convert_names(conversion_file: str, 
+                  dataset: xr.Dataset,
+                  data_variable: str,
+                  convert_seasons_and_terms: bool = False):
+
+    # Load conversion dictionaries
+    with open(conversion_file, 'rb') as f:
+        converter = pickle.load(f)
+
+    # Convert
+    new_dataset = convert_coordname_elements(dataset, 
+                                             converter['coord_names'], 
+                                             converter['coord_element_names'],
+                                            False)   
+
+    # Convert weeks and hours
+    if convert_seasons_and_terms:
+        new_dataset = (
+            new_dataset
+            .assign_coords(S=converter['week_to_seasons'])
+            .assign_coords(T=converter['hour_to_terms'])
+        )
+
+    # Test that we did not mess something up
+    before = np.nan_to_num(dataset.get(data_variable).data)
+    after = np.nan_to_num(new_dataset.get(data_variable).data)
+    assert np.all(after == before), 'Values are not the same after conversion!'
+
+    return dataset, new_dataset
+
+def transform_xrdata(xarray: xr.Dataset,
+                  data_variable: str,
+                  selection: Tuple[dict, None] = None,
+                  sets_to_sum: Tuple[str, list, None] = None):
+
+    output = xarray
+    if sets_to_sum != None:
+        output = (
+            output
+            .get(data_variable)
+            .sum(sets_to_sum)
+        )
+
+    if selection != None:
+        output = output.sel(selection)
+
+    return output

RawDataProcessing/Modules/Submodules/aauvarmeplan2021.py → RawDataProcessing/Modules/aauvarmeplan2021.py

@@ -47,7 +47,7 @@
 import numpy as np
 import geopandas as gpd
 import os
-from Modules.Submodules.municipal_template import DataContainer
+from Submodules.municipal_template import DataContainer
 
 style = 'report'
 
@@ -127,7 +127,7 @@ def __init__(self, scenario: str = 'SUM_GWh_uden_besp') -> None:
 
             ## Sum to scenario
             f = f.pivot_table(index=['year', 'municipality', 'user'], values=scenario, 
-                              aggfunc=lambda x: np.sum(x)*1e3) # To MWh
+                              aggfunc=lambda x: np.sum(x)*1e3, observed=False) # To MWh
 
             self.DH = pd.concat((self.DH, f))
 
@@ -146,26 +146,36 @@ def __init__(self, scenario: str = 'SUM_GWh_uden_besp') -> None:
         self.IND = temp.to_xarray()
 
 
+def main():
+    # 1.1 Get formatted Varmeplan2021 Data
+    data_format = VPDK21()
+    data_format.DH.to_netcdf('Data/AAU Kommuneplan/districtheat_exo_heatdem.nc')
+    data_format.IND.to_netcdf('Data/AAU Kommuneplan/industry_exo_heatdem.nc')
+
 if __name__ == '__main__':
-    data = DataContainer()
-    VP = VPDK21()
-    data.muni = data.muni.merge(VP.DH)
-    data.muni = data.muni.merge(VP.IND)
+    main()    
+
+    merge_example = False
+    if merge_example:
+        data = DataContainer()
+        VP = VPDK21()
+        data.muni = data.muni.merge(VP.DH)
+        data.muni = data.muni.merge(VP.IND)
 
-    # Plotting geometry
-    geo = data.get_polygons('muni')
-    for user in data.muni.user.data:
-        if user == 'district_heating' or user == 'individual':
-            geo[user] = data.muni.heat_demand_mwh.sel(user=user,
-                                                    year=2019).data
-            unit = 'MWh'
-        else:
-            geo[user] = data.muni.heat_demand_normalised.sel(user=user,
-                                                    year=2019).data
-            unit = 'normalised'
-
-        geo.plot(
-                    column=user,
-                    cmap='coolwarm',  # Add the cmap parameter here
-                    legend=True       # Add the legend parameter here
-                ).set_title(user + f' Heat Demand ({unit})')
+        # Plotting geometry
+        geo = data.get_polygons('muni')
+        for user in data.muni.user.data:
+            if user == 'district_heating' or user == 'individual':
+                geo[user] = data.muni.heat_demand_mwh.sel(user=user,
+                                                        year=2019).data
+                unit = 'MWh'
+            else:
+                geo[user] = data.muni.heat_demand_normalised.sel(user=user,
+                                                        year=2019).data
+                unit = 'normalised'
+            
+            geo.plot(
+                        column=user,
+                        cmap='coolwarm',  # Add the cmap parameter here
+                        legend=True       # Add the legend parameter here
+                    ).set_title(user + f' Heat Demand ({unit})')

RawDataProcessing/Modules/createDE.py

@@ -21,6 +21,7 @@
 import numpy as np
 import click
 import pickle
+from Submodules.utils import convert_names, transform_xrdata 
 import xarray as xr
 from typing import Tuple
 from pybalmorel import IncFile
@@ -51,61 +52,6 @@
 ###      1. Make .inc Files         ###
 ### ------------------------------- ###
 
-def convert_names(conversion_file: str, el_dataset: str):
-    """_summary_
-
-    Args:
-        conversion_file (str): _description_
-        el_dataset (str): _description_
-
-    Returns:
-        _type_: _description_
-    """
-    # Load dataset
-    dataset = xr.load_dataset(el_dataset)
-
-    # Load conversion dictionaries
-    with open(conversion_file, 'rb') as f:
-        converter = pickle.load(f)
-
-    # Convert
-    new_dataset = convert_coordname_elements(dataset, 
-                                             converter['coord_names'], 
-                                             converter['coord_element_names'],
-                                            False)   
-
-    # Convert weeks and hours
-    new_dataset = (
-        new_dataset
-        .assign_coords(S=converter['week_to_seasons'])
-        .assign_coords(T=converter['hour_to_terms'])
-    )
-
-    # Test that we did not mess something up
-    before = np.nan_to_num(dataset.electricity_demand_mwh.data)
-    after = np.nan_to_num(new_dataset.electricity_demand_mwh.data)
-    assert np.all(after == before), 'Values are not the same after conversion!'
-
-    return dataset, new_dataset
-
-def transform_xrdata(xarray: xr.Dataset,
-                  data_variable: str,
-                  selection: Tuple[dict, None] = None,
-                  sets_to_sum: Tuple[str, list, None] = None):
-
-    output = xarray
-    if sets_to_sum != None:
-        output = (
-            output
-            .get(data_variable)
-            .sum(sets_to_sum)
-        )
-
-    if selection != None:
-        output = output.sel(selection)
-
-    return output
-
 # Main function
 @click.command()
 @click.option("--conversion-file", type=str, required=True, help="The conversion dictionary")
@@ -126,8 +72,11 @@ def main(conversion_file: str,
         None
     """
 
-    # 1.1 Format Dataset
-    dataset, new_dataset = convert_names(conversion_file, el_dataset)
+    # 1.1 Format Dataset    
+    ## Load dataset
+    el_dataset = xr.load_dataset(el_dataset) # converts from string path to xr.Dataset
+    dataset, new_dataset = convert_names(conversion_file, el_dataset, 
+                                         'electricity_demand_mwh', convert_seasons_and_terms=True)
 
     if show_difference:
         print('Before: \n', dataset, '\n\n')
@@ -143,6 +92,7 @@ def main(conversion_file: str,
         suffix=textwrap.dedent(""";
         DE(YYY,RRR,DEUSER)=DE1(RRR,DEUSER,YYY);
         DE('2050',RRR,DEUSER) = DE('2023', RRR, DEUSER);
+        DE1(RRR,DEUSER,YYY) = 0;
         """),
         path=out_path
     )