Add Validation Charts

workflow/rules/postprocess.smk

@@ -71,4 +71,17 @@ rule summarise_results:
     log:
         log = 'results/{scenario}/logs/summarise_results.log'
     shell: 
-        'python workflow/scripts/osemosys_global/summarise_results.py 2> {log}'
+        'python workflow/scripts/osemosys_global/summarise_results.py 2> {log}'
+
+rule calculate_carbon_intensity:
+    message:
+        "Calculating Carbon Intensity..."
+    input:
+        production_by_technology = "results/{scenario}/results/ProductionByTechnologyAnnual.csv",
+        annual_emissions = "results/{scenario}/results/AnnualEmissions.csv",
+    output:
+        emission_intensity = "results/{scenario}/results/AnnualEmissionIntensity.csv",
+    log:
+        log = 'results/{scenario}/logs/calcualte_carbon_intensity.log'
+    script: 
+        "../scripts/osemosys_global/carbon_intensity.py"

workflow/rules/validate.smk

@@ -0,0 +1,119 @@
+"""Rules for validating historical results"""
+
+###
+# constants
+###
+
+# used for target rules
+CAPACITY_VALIDATION = ["ember", "irena", "eia"]
+GENERATION_VALIDATION = ["ember", "irena", "eia"]
+EMISSION_VALIDATION = ["ember", "climatewatch"]
+EMISSION_INTENSITY_VALIDATION = ["ember"]
+
+###
+# capacity
+###
+
+def capacity_validation_data(wildcards):
+    if wildcards.datasource == "ember":
+        return "resources/data/ember_yearly_electricity_data.csv"
+    elif wildcards.datasource == "irena":
+        return "resources/data/validation/irena_capacity.csv"
+    elif wildcards.datasource == "eia":
+        return "resources/data/validation/eia_capacity.json"
+
+def capacity_options(wildcards):
+    if wildcards.datasource == "irena":
+        return {"iso_codes": "resources/data/validation/iso.csv"}
+    else:
+        return []
+
+rule validate_capacity:
+    message: "Validating capacity against {wildcards.datasource}"
+    input:
+        unpack(capacity_options),
+        validation_data = capacity_validation_data,
+        og_result = "results/{scenario}/results/TotalCapacityAnnual.csv",
+    params:
+        result_dir="results/{scenario}/results",
+        variable="capacity"
+    output:
+        expand("results/{{scenario}}/validation/{country}/capacity/{{datasource}}.png", country=COUNTRIES)
+    script:
+        "../scripts/osemosys_global/validation/main.py"
+
+###
+# generation
+###
+
+def generation_validation_data(wildcards):
+    if wildcards.datasource == "ember":
+        return "resources/data/ember_yearly_electricity_data.csv"
+    elif wildcards.datasource == "irena":
+        return "resources/data/validation/irena_generation.csv"
+    elif wildcards.datasource == "eia":
+        return "resources/data/validation/eia_generation.json"
+
+def generation_options(wildcards):
+    if wildcards.datasource == "irena":
+        return {"iso_codes": "resources/data/validation/iso.csv"}
+    else:
+        return []
+
+rule validate_generation:
+    message: "Validating generation against {wildcards.datasource}"
+    input:
+        unpack(generation_options),
+        validation_data = generation_validation_data,
+        og_result = "results/{scenario}/results/ProductionByTechnologyAnnual.csv",
+    params:
+        result_dir="results/{scenario}/results",
+        variable="generation"
+    output:
+        expand("results/{{scenario}}/validation/{country}/generation/{{datasource}}.png", country=COUNTRIES)
+    script:
+        "../scripts/osemosys_global/validation/main.py"
+
+###
+# emissions
+###
+
+def emission_validation_data(wildcards):
+    if wildcards.datasource == "ember":
+        return "resources/data/ember_yearly_electricity_data.csv"
+    elif wildcards.datasource == "climatewatch":
+        return "resources/data/validation/climatewatch.csv"
+
+rule validate_emissions:
+    message: "Validating emissions against {wildcards.datasource}"
+    input:
+        validation_data = emission_validation_data,
+        og_result = "results/{scenario}/results/AnnualEmissions.csv"
+    params:
+        result_dir="results/{scenario}/results",
+        variable="emissions"
+    output:
+        expand("results/{{scenario}}/validation/{country}/emissions/{{datasource}}.png", country=COUNTRIES)
+    script:
+        "../scripts/osemosys_global/validation/main.py"
+
+###
+# emission intensity
+###
+
+def emission_intensity_validation_data(wildcards):
+    if wildcards.datasource == "ember":
+        return "resources/data/ember_yearly_electricity_data.csv"
+
+rule validate_emission_intensity:
+    message: "Validating emission intensity against {wildcards.datasource}"
+    input:
+        validation_data = emission_intensity_validation_data,
+        og_result = "results/{scenario}/results/AnnualEmissionIntensity.csv"
+    params:
+        result_dir="results/{scenario}/results",
+        variable="emission_intensity"
+    output:
+        expand("results/{{scenario}}/validation/{country}/emission_intensity/{{datasource}}.png", country=COUNTRIES)
+    script:
+        "../scripts/osemosys_global/validation/main.py"

workflow/scripts/osemosys_global/carbon_intensity.py

@@ -0,0 +1,74 @@
+"""Calcualtes carbon intensity
+
+Since emissions are reported at country levels, carbon intensities are also 
+reported at country levels. 
+"""
+
+import pandas as pd
+
+
+def format_production(production_by_technology_annual: pd.DataFrame) -> pd.DataFrame:
+
+    df = production_by_technology_annual.copy()
+
+    df = df[df.TECHNOLOGY.str.startswith("PWR")]
+    df["TECH"] = df.TECHNOLOGY.str[3:6]
+    df["CTRY"] = df.TECHNOLOGY.str[6:9]
+    df = df[~(df.TECH.isin(["BAT", "TRN"]))]
+    df = df.drop(columns=["TECHNOLOGY", "FUEL", "TECH"])
+    return df.groupby(["REGION", "CTRY", "YEAR"]).sum()
+
+
+def format_emissions(annual_emissions: pd.DataFrame) -> pd.DataFrame:
+
+    df = annual_emissions.copy()
+
+    df["CTRY"] = df.EMISSION.str[3:6]
+    return df.groupby(["REGION", "EMISSION", "CTRY", "YEAR"]).sum()
+
+
+def calculate_emission_intensity(
+    production: pd.DataFrame, emissions: pd.DataFrame
+) -> pd.DataFrame:
+
+    p = production.copy().rename(columns={"VALUE": "PRODUCTION_PJ"})
+    e = emissions.copy().rename(columns={"VALUE": "EMISSIONS_MT"})
+
+    df = p.join(e).fillna(0)
+    df = df.droplevel("CTRY")  # emission retains country
+
+    # 1 PJ (1000 TJ / PJ) (1000 GJ / TJ) (1000 MJ / GJ) (1000 kJ / MJ) (hr / 3600sec)
+    df["PRODUCTION_kwh"] = (
+        df.PRODUCTION_PJ.mul(1000).mul(1000).mul(1000).mul(1000).div(3600)
+    )
+
+    # 1 MT (1,000,000 T / MT) (1000 kg / T) (1000g / kg)
+    df["EMISSIONS_g"] = df.EMISSIONS_MT.mul(1000000).mul(1000).mul(1000)
+
+    # intensity in g/kwh
+    df["VALUE"] = df.EMISSIONS_g.div(df.PRODUCTION_kwh)
+
+    return df["VALUE"].to_frame()
+
+
+if __name__ == "__main__":
+    if "snakemake" in globals():
+        production_by_technology_annual_csv = snakemake.input.production_by_technology
+        annual_emissions_csv = snakemake.input.annual_emissions
+        save = snakemake.output.emission_intensity
+    else:
+        production_by_technology_annual_csv = (
+            "results/India/results/ProductionByTechnologyAnnual.csv"
+        )
+        annual_emissions_csv = "results/India/results/AnnualEmissions.csv"
+        save = "results/India/results/AnnualEmissionIntensity.csv"
+
+    production_by_technology_annual = pd.read_csv(production_by_technology_annual_csv)
+    annual_emissions = pd.read_csv(annual_emissions_csv)
+
+    prod = format_production(production_by_technology_annual)
+    emissions = format_emissions(annual_emissions)
+
+    emission_intensity = calculate_emission_intensity(prod, emissions).round(2)
+
+    emission_intensity.to_csv(save, index=True)

workflow/scripts/osemosys_global/validation/climate_watch.py

@@ -0,0 +1,40 @@
+"""Data handeling for Climate Watch validation
+
+https://www.climatewatchdata.org/ghg-emissions
+"""
+
+import pandas as pd
+
+###
+# public functions
+###
+
+
+def get_cw_emissions(csv_file: str, **kwargs) -> pd.DataFrame:
+    """Gets Climate Watch emissions data"""
+    df = _read_cw_data(csv_file)
+    return _format_cw_data(df)
+
+
+###
+# private functions
+###
+
+
+def _read_cw_data(csv_file: str) -> pd.DataFrame:
+    """Reads climate watch data
+
+    https://www.climatewatchdata.org/ghg-emissions?end_year=2021&gases=all-ghg&sectors=electricity-heat&start_year=1990
+    """
+    return pd.read_csv(csv_file, skipfooter=2, engine="python")
+
+
+def _format_cw_data(cw: pd.DataFrame) -> pd.DataFrame:
+    df = cw.copy()
+    df = df.drop(columns=["Country/Region", "unit"]).rename(columns={"iso": "EMISSION"})
+    df = df.melt(id_vars=["EMISSION"], var_name="YEAR", value_name="VALUE")
+    df = df.fillna(0).replace("false", 0)
+    df["YEAR"] = df.YEAR.astype(int)
+    df["VALUE"] = df.VALUE.astype(float)
+    df["REGION"] = "GLOBAL"
+    return df.set_index(["REGION", "EMISSION", "YEAR"])