notebooks/resources.py.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0b71c0bc-9fb6-47b9-b75e-ea6d6df7fd83",
   "metadata": {},
   "outputs": [],
   "source": [
    "import pypsa\n",
    "import yaml\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import geopandas as gpd\n",
    "import xarray as xr\n",
    "import cartopy.crs as ccrs\n",
    "import cartopy\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "from matplotlib.lines import Line2D\n",
    "import matplotlib.colors as mcolors\n",
    "from matplotlib.patches import Circle, Patch\n",
    "from matplotlib.legend_handler import HandlerPatch\n",
    "\n",
    "from pypsa.descriptors import get_switchable_as_dense as as_dense\n",
    "\n",
    "from shapely import wkt\n",
    "import sys, os\n",
    "\n",
    "\n",
    "PATH = \"../workflows/\"\n",
    "SCRIPTS_PATH = \"pypsa-eur-sec/scripts/\"\n",
    "\n",
    "sys.path.append(os.path.join(PATH, SCRIPTS_PATH))\n",
    "from plot_summary import rename_techs\n",
    "from plot_network import assign_location\n",
    "from helper import override_component_attrs\n",
    "from build_gas_input_locations import build_gas_input_locations, load_bus_regions\n",
    "\n",
    "plt.style.use([\"bmh\", \"matplotlibrc\"])\n",
    "xr.set_options(display_style=\"html\")\n",
    "\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ce4faa73-c438-4193-8274-2046d9d2d812",
   "metadata": {},
   "outputs": [],
   "source": [
    "OUTPUT = \"../results/graphics-20221227\"\n",
    "\n",
    "CLUSTERS = 181\n",
    "LL = \"opt\"\n",
    "OPTS = \"Co2L0-3H-T-H-B-I-A-solar+p3-linemaxext10\"\n",
    "RUN = \"20221227-gas\"\n",
    "SCENARIO = f\"elec_s_{CLUSTERS}_lv{LL}__{OPTS}_2050\"\n",
    "OVERRIDES = PATH + \"pypsa-eur-sec/data/override_component_attrs\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "618d5fc8-ef0b-4980-87a0-7dbae216de05",
   "metadata": {},
   "outputs": [],
   "source": [
    "with open(PATH + \"pypsa-eur-sec/config.main.yaml\") as file:\n",
    "    config = yaml.safe_load(file)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "28fcd3d0-a18d-4dae-83dd-7b2883117a91",
   "metadata": {},
   "outputs": [],
   "source": [
    "fn = f\"{PATH}/pypsa-eur/resources/regions_onshore_elec_s_{CLUSTERS}.geojson\"\n",
    "nodes = gpd.read_file(fn).set_index(\"name\")\n",
    "\n",
    "fn = f\"{PATH}/pypsa-eur/resources/regions_offshore_elec_s_{CLUSTERS}.geojson\"\n",
    "offnodes = gpd.read_file(fn).set_index(\"name\")\n",
    "\n",
    "fn = f\"{PATH}/pypsa-eur/resources/country_shapes.geojson\"\n",
    "cts = gpd.read_file(fn).set_index(\"name\")\n",
    "\n",
    "regions = pd.concat(\n",
    "    [\n",
    "        gpd.read_file(f\"{PATH}/pypsa-eur/resources/regions_onshore.geojson\"),\n",
    "        gpd.read_file(f\"{PATH}/pypsa-eur/resources/regions_offshore.geojson\"),\n",
    "    ]\n",
    ")\n",
    "regions = regions.dissolve(\"name\")\n",
    "\n",
    "fn = f\"{PATH}/pypsa-eur/resources/regions_onshore.geojson\"\n",
    "onregions = gpd.read_file(fn).set_index(\"name\")\n",
    "\n",
    "fn = f\"{PATH}/pypsa-eur/resources/regions_onshore.geojson\"\n",
    "offregions = gpd.read_file(fn).set_index(\"name\")\n",
    "\n",
    "epsg = 3035\n",
    "regions[\"Area\"] = regions.to_crs(epsg=epsg).area.div(1e6)\n",
    "onregions[\"Area\"] = onregions.to_crs(epsg=epsg).area.div(1e6)\n",
    "offregions[\"Area\"] = offregions.to_crs(epsg=epsg).area.div(1e6)\n",
    "nodes[\"Area\"] = nodes.to_crs(epsg=epsg).area.div(1e6)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9df4a814-45cf-4706-b4ac-d4056ecfcc25",
   "metadata": {},
   "outputs": [],
   "source": [
    "europe_shape = nodes.dissolve()\n",
    "europe_shape.index = [\"EU\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f0635d35-68f1-47fd-86e6-7bcc10eb25c4",
   "metadata": {},
   "outputs": [],
   "source": [
    "minx, miny, maxx, maxy = europe_shape.explode(ignore_index=True).total_bounds\n",
    "BOUNDARIES = [minx, maxx - 4, miny, maxy]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8570150a-ea32-4e8a-9599-3b05c28cf935",
   "metadata": {},
   "outputs": [],
   "source": [
    "overrides = override_component_attrs(OVERRIDES)\n",
    "fn = f\"{PATH}/pypsa-eur-sec/results/{RUN}/prenetworks/{SCENARIO}.nc\"\n",
    "n = pypsa.Network(fn, override_component_attrs=overrides)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "49dcc78f-c79c-4f79-a658-295189336a43",
   "metadata": {},
   "outputs": [],
   "source": [
    "unique_link_carriers = n.links.carrier.unique()\n",
    "GAS_NETWORK = \"gas pipeline\" in unique_link_carriers\n",
    "H2_NETWORK = any(\"H2 pipeline\" in ulc for ulc in unique_link_carriers)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "51b7a780-c45f-4a5c-88d2-ad5cfd06a62b",
   "metadata": {},
   "source": [
    "## Utilities"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cb272cdb-ca36-49a1-915f-acb92106f19c",
   "metadata": {},
   "outputs": [],
   "source": [
    "def rename_techs_tyndp(tech):\n",
    "    tech = rename_techs(tech)\n",
    "    # if \"heat pump\" in tech or \"resistive heater\" in tech:\n",
    "    #    return \"power-to-heat\"\n",
    "    # elif tech in [\"H2 Electrolysis\", \"methanation\", \"helmeth\", \"H2 liquefaction\"]:\n",
    "    #    return \"power-to-gas\"\n",
    "    if tech == \"H2\":\n",
    "        return \"H2 storage\"\n",
    "    # elif tech in [\"OCGT\", \"CHP\", \"gas boiler\", \"H2 Fuel Cell\"]:\n",
    "    #    return \"gas-to-power/heat\"\n",
    "    # elif \"solar\" in tech:\n",
    "    #    return \"solar\"\n",
    "    elif tech == \"Fischer-Tropsch\":\n",
    "        return \"power-to-liquid\"\n",
    "    elif \"offshore wind\" in tech:\n",
    "        return \"offshore wind\"\n",
    "    #    if \"heat pump\" in tech:\n",
    "    #        return \"heat pump\"\n",
    "    elif tech == \"gas\":\n",
    "        return \"fossil gas\"\n",
    "    # elif \"CC\" in tech or \"sequestration\" in tech:\n",
    "    #    return \"CCS\"\n",
    "    elif tech in [\"industry electricity\", \"agriculture electricity\"]:\n",
    "        return \"industry electricity\"\n",
    "    elif \"oil emissions\" in tech:\n",
    "        return \"oil emissions\"\n",
    "    else:\n",
    "        return tech"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6e1923ad-6d04-4bf6-a14d-c9bcf4b1e341",
   "metadata": {},
   "outputs": [],
   "source": [
    "class HandlerCircle(HandlerPatch):\n",
    "    \"\"\"\n",
    "    Legend Handler used to create circles for legend entries.\n",
    "\n",
    "    This handler resizes the circles in order to match the same dimensional\n",
    "    scaling as in the applied axis.\n",
    "    \"\"\"\n",
    "\n",
    "    def create_artists(\n",
    "        self, legend, orig_handle, xdescent, ydescent, width, height, fontsize, trans\n",
    "    ):\n",
    "        fig = legend.get_figure()\n",
    "        ax = legend.axes\n",
    "\n",
    "        unit = np.diff(ax.transData.transform([(0, 0), (1, 1)]), axis=0)[0][1]\n",
    "        radius = orig_handle.get_radius() * unit * (72 / fig.dpi)\n",
    "        center = 5 - xdescent, 3 - ydescent\n",
    "        p = plt.Circle(center, radius)\n",
    "        self.update_prop(p, orig_handle, legend)\n",
    "        p.set_transform(trans)\n",
    "        return [p]\n",
    "\n",
    "\n",
    "def add_legend_circles(\n",
    "    ax, sizes, labels, scale=1, srid=None, patch_kw={}, legend_kw={}\n",
    "):\n",
    "    if srid is not None:\n",
    "        area_correction = projected_area_factor(ax, n.srid) ** 2\n",
    "        sizes = [s * area_correction for s in sizes]\n",
    "\n",
    "    handles = make_legend_circles_for(sizes, scale, **patch_kw)\n",
    "\n",
    "    legend = ax.legend(\n",
    "        handles, labels, handler_map={Circle: HandlerCircle()}, **legend_kw\n",
    "    )\n",
    "\n",
    "    ax.add_artist(legend)\n",
    "\n",
    "\n",
    "def add_legend_lines(ax, sizes, labels, scale=1, patch_kw={}, legend_kw={}):\n",
    "    handles = [Line2D([0], [0], linewidth=s / scale, **patch_kw) for s in sizes]\n",
    "\n",
    "    legend = ax.legend(handles, labels, **legend_kw)\n",
    "\n",
    "    ax.add_artist(legend)\n",
    "\n",
    "\n",
    "def add_legend_patch(ax, colors, labels, patch_kw={}, legend_kw={}):\n",
    "    handles = [Patch(facecolor=c, **patch_kw) for c in colors]\n",
    "\n",
    "    legend = ax.legend(handles, labels, **legend_kw)\n",
    "\n",
    "    ax.add_artist(legend)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "16d92671-05d5-499b-b9dd-78c6e4c790a6",
   "metadata": {},
   "source": [
    "## Original Gas and Electricity Network"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "09d0e095-b9ad-4d27-92dd-07611d9ee71e",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_original_gas_network(network):\n",
    "    n = network.copy()\n",
    "\n",
    "    lw_factor = 20e3\n",
    "    color = \"indianred\"\n",
    "\n",
    "    n.mremove(\"Bus\", n.buses.index[n.buses.carrier != \"AC\"])\n",
    "    n.mremove(\"Link\", n.links.index[~n.links.carrier.str.contains(\"gas pipeline\")])\n",
    "\n",
    "    n.links.bus0 = n.links.bus0.str.replace(\" gas\", \"\")\n",
    "    n.links.bus1 = n.links.bus1.str.replace(\" gas\", \"\")\n",
    "\n",
    "    crs = ccrs.EqualEarth()\n",
    "\n",
    "    fig, ax = plt.subplots(figsize=(7, 6), subplot_kw={\"projection\": crs})\n",
    "\n",
    "    n.plot(\n",
    "        geomap=True,\n",
    "        ax=ax,\n",
    "        bus_sizes=0.005,\n",
    "        bus_colors=\"k\",\n",
    "        link_colors=color,\n",
    "        link_widths=n.links.p_nom / lw_factor,\n",
    "        branch_components=[\"Link\"],\n",
    "    )\n",
    "\n",
    "    handles = []\n",
    "    labels = []\n",
    "\n",
    "    for s in (100, 50, 20):\n",
    "        handles.append(plt.Line2D([0], [0], color=color, linewidth=s * 1e3 / lw_factor))\n",
    "        labels.append(\"{} GW\".format(s))\n",
    "    l1_1 = ax.legend(\n",
    "        handles,\n",
    "        labels,\n",
    "        loc=\"upper left\",\n",
    "        bbox_to_anchor=(0.05, 1.01),\n",
    "        frameon=False,\n",
    "        labelspacing=0.8,\n",
    "        handletextpad=1.5,\n",
    "        title=\"Today's gas transmission\",\n",
    "    )\n",
    "    ax.add_artist(l1_1)\n",
    "\n",
    "    fig.savefig(f\"{OUTPUT}/gas-network-today-map.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8975454b-4200-42bb-acbe-0ebbabd9b3fd",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_original_gas_network(n)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "98f32688-61c0-4cf4-b56b-8ae7007c0f15",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_original_electricity_network(network):\n",
    "    n = network.copy()\n",
    "\n",
    "    lw_factor = 3e3\n",
    "\n",
    "    n.mremove(\"Bus\", n.buses.index[n.buses.carrier != \"AC\"])\n",
    "    n.mremove(\"Link\", n.links.index[~n.links.carrier.str.contains(\"DC\")])\n",
    "\n",
    "    crs = ccrs.EqualEarth()\n",
    "\n",
    "    fig, ax = plt.subplots(figsize=(7, 6), subplot_kw={\"projection\": crs})\n",
    "\n",
    "    n.plot(\n",
    "        geomap=True,\n",
    "        ax=ax,\n",
    "        bus_sizes=0.005,\n",
    "        bus_colors=\"k\",\n",
    "        link_colors=n.links.p_nom.apply(\n",
    "            lambda x: \"darkseagreen\" if x > 0 else \"skyblue\"\n",
    "        ),\n",
    "        link_widths=1.5,\n",
    "        line_widths=n.lines.s_nom / lw_factor,\n",
    "    )\n",
    "\n",
    "    handles = []\n",
    "    labels = []\n",
    "\n",
    "    for s in (20, 10):\n",
    "        handles.append(\n",
    "            plt.Line2D([0], [0], color=\"rosybrown\", linewidth=s * 1e3 / lw_factor)\n",
    "        )\n",
    "        labels.append(\"{} GW HVAC today\".format(s))\n",
    "    l1_1 = ax.legend(\n",
    "        handles,\n",
    "        labels,\n",
    "        loc=\"upper left\",\n",
    "        bbox_to_anchor=(0.025, 1.01),\n",
    "        frameon=False,\n",
    "        labelspacing=0.8,\n",
    "        handletextpad=1.5,\n",
    "        # title='Today\\'s AC transmission'\n",
    "    )\n",
    "    ax.add_artist(l1_1)\n",
    "\n",
    "    handles = [\n",
    "        Line2D([0], [0], color=\"darkseagreen\", lw=1.5),\n",
    "        Line2D([0], [0], color=\"skyblue\", lw=1.5),\n",
    "    ]\n",
    "    l2 = ax.legend(\n",
    "        handles,\n",
    "        [\"HVDC existing\", \"HVDC planned\"],\n",
    "        frameon=False,\n",
    "        loc=[0.04, 0.78],\n",
    "        labelspacing=0.8,\n",
    "        handletextpad=1.5,\n",
    "    )\n",
    "\n",
    "    ax.add_artist(l2)\n",
    "\n",
    "    fig.savefig(f\"{OUTPUT}/electricity-network-today-map.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "652fd996-9b30-44c6-be2a-160d4bc49bf3",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_original_electricity_network(n)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "df00bf4e-6585-408b-81db-3b7c5a81f0f8",
   "metadata": {},
   "source": [
    "## Stats: Per-Country Potentials"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "07a8cde8-7e99-494e-b8e7-a23c70f1dc2c",
   "metadata": {},
   "outputs": [],
   "source": [
    "def country_potentials(n):\n",
    "    return (\n",
    "        n.generators.p_nom_max.groupby([n.generators.carrier, n.generators.bus.str[:2]])\n",
    "        .sum()\n",
    "        .div(1e3)\n",
    "        .unstack(0)\n",
    "        .round(1)\n",
    "        .replace([np.inf, -np.inf], np.nan)\n",
    "        .dropna(axis=1, how=\"all\")\n",
    "        .fillna(0.0)\n",
    "    )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a06a18d9-3b50-4bf8-b87c-4933b4fa715a",
   "metadata": {},
   "outputs": [],
   "source": [
    "p = country_potentials(n)\n",
    "p"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9a1c519b-8ecf-49bb-906c-3d7ed3a684b2",
   "metadata": {},
   "outputs": [],
   "source": [
    "p.to_csv(f\"../results/tables/country_potentials.csv\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5b57db3e-a2ef-447c-a8a5-96d8f2c323df",
   "metadata": {},
   "source": [
    "## Maps: regional demands"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fb475e16-2bcc-4473-acf6-ba01deb91988",
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_oil_demand(with_methanol=False):\n",
    "    fn = f\"{PATH}/pypsa-eur-sec/resources/industrial_energy_demand_elec_s_{CLUSTERS}_2050.csv\"\n",
    "    industrial_demand = pd.read_csv(fn, index_col=0)\n",
    "\n",
    "    fn = f\"{PATH}/pypsa-eur-sec/resources/pop_layout_elec_s_{CLUSTERS}.csv\"\n",
    "    pop_layout = pd.read_csv(fn, index_col=0)\n",
    "\n",
    "    fn = f\"{PATH}/pypsa-eur-sec/resources/pop_weighted_energy_totals_s_{CLUSTERS}.csv\"\n",
    "    nodal_energy_totals = pd.read_csv(fn, index_col=0)\n",
    "\n",
    "    oil = [\n",
    "        \"total international aviation\",\n",
    "        \"total domestic aviation\",\n",
    "        \"total agriculture machinery\",\n",
    "    ]\n",
    "    if with_methanol:\n",
    "        oil.append(\"total domestic navigation\")\n",
    "\n",
    "    demand = industrial_demand[\"naphtha\"] + nodal_energy_totals[oil].sum(axis=1)\n",
    "\n",
    "    if with_methanol:\n",
    "        fn = f\"{PATH}/pypsa-eur-sec/resources/shipping_demand_s_{CLUSTERS}.csv\"\n",
    "        efficiency = (\n",
    "            config[\"sector\"][\"shipping_oil_efficiency\"]\n",
    "            / config[\"sector\"][\"shipping_methanol_efficiency\"]\n",
    "        )\n",
    "        demand += pd.read_csv(fn, index_col=0).squeeze() * efficiency\n",
    "\n",
    "    return demand"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3b568138-39b0-44b8-84e1-7d98a67a5234",
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_process_emission():\n",
    "    fn = f\"{PATH}/pypsa-eur-sec/resources/industrial_energy_demand_elec_s_{CLUSTERS}_2050.csv\"\n",
    "    industrial_demand = pd.read_csv(fn, index_col=0)\n",
    "\n",
    "    return industrial_demand[\"process emission\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d23d25a5-7763-40a4-b9ef-13dabc1080ad",
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_biomass_demand():\n",
    "    fn = f\"{PATH}/pypsa-eur-sec/resources/industrial_energy_demand_elec_s_{CLUSTERS}_2050.csv\"\n",
    "    industrial_demand = pd.read_csv(fn, index_col=0)\n",
    "\n",
    "    return industrial_demand[\"solid biomass\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5ffe405c-3d79-4ac7-94bd-5967e1b31005",
   "metadata": {},
   "outputs": [],
   "source": [
    "title = {\n",
    "    \"electricity\": \"Electricity Demand [TWh/a]\",\n",
    "    \"H2\": \"Hydrogen Demand [TWh/a]\",\n",
    "    \"heat\": \"Heat Demand [TWh/a]\",\n",
    "    \"solid biomass\": \"Solid Biomass Demand [TWh/a]\",\n",
    "    \"gas\": \"Methane Demand [TWh/a]\",\n",
    "    \"oil\": \"Liquid Hydrocarbons Demand [TWh/a]\",\n",
    "    \"process emission\": \"Process Emissions [MtCO2/a]\",\n",
    "}\n",
    "\n",
    "cmap = {\n",
    "    \"electricity\": \"Blues\",\n",
    "    \"H2\": \"RdPu\",\n",
    "    \"heat\": \"Reds\",\n",
    "    \"solid biomass\": \"Greens\",\n",
    "    \"gas\": \"Oranges\",\n",
    "    \"oil\": \"Greys\",\n",
    "    \"process emission\": \"Greys\",\n",
    "}\n",
    "\n",
    "regex = {\n",
    "    \"electricity\": r\"(electricity|EV)\",\n",
    "    \"H2\": r\"(H2|fuel cell)\",\n",
    "    \"heat\": r\"heat\",\n",
    "    \"solid biomass\": r\"biomass\",\n",
    "    \"oil\": r\"( oil|naphtha|kerosene|methanol)\",\n",
    "    \"gas\": r\"gas\",\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ca066708-36f1-4b86-8f3c-a5b2592e9fd1",
   "metadata": {},
   "outputs": [],
   "source": [
    "demand = as_dense(n, \"Load\", \"p_set\").div(1e6)  # TWh\n",
    "demand_grouped = demand.groupby(\n",
    "    [n.loads.carrier, n.loads.bus.map(n.buses.location)], axis=1\n",
    ").sum()\n",
    "demand_by_region = (n.snapshot_weightings.generators @ demand_grouped).unstack(level=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5567fdaa-236c-4647-9e1a-9b38af756ed6",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_regional_demands(df, geodf, carrier, series=None, vmax=None, vmin=0):\n",
    "    if series is None:\n",
    "        series = df.filter(regex=regex[carrier]).sum(axis=1).drop(\"EU\")\n",
    "\n",
    "    proj = ccrs.EqualEarth()\n",
    "    geodf = geodf.to_crs(proj.proj4_init)\n",
    "\n",
    "    fig, ax = plt.subplots(figsize=(7, 7), subplot_kw={\"projection\": proj})\n",
    "\n",
    "    geodf.plot(\n",
    "        ax=ax,\n",
    "        column=series,\n",
    "        # transform=ccrs.PlateCarree(),\n",
    "        cmap=cmap[carrier],\n",
    "        linewidths=0,\n",
    "        legend=True,\n",
    "        vmin=vmin,\n",
    "        vmax=vmax,\n",
    "        legend_kwds={\"label\": title[carrier], \"shrink\": 0.7, \"extend\": \"max\"},\n",
    "    )\n",
    "\n",
    "    ax.add_feature(cartopy.feature.COASTLINE.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "    ax.add_feature(cartopy.feature.BORDERS.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "\n",
    "    ax.set_frame_on(False)\n",
    "    ax.set_facecolor(\"white\")\n",
    "\n",
    "    plt.savefig(f\"{OUTPUT}/demand-map-{carrier}.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2c305aef-ff77-486a-a0cc-a72b834da221",
   "metadata": {},
   "outputs": [],
   "source": [
    "report = {}\n",
    "report[\"demand\"] = pd.concat(\n",
    "    {\n",
    "        c: demand_by_region.filter(regex=regex[c]).sum(axis=1).drop(\"EU\")\n",
    "        for c in [\"electricity\", \"H2\", \"heat\"]\n",
    "    }\n",
    "    | {\"oil\": get_oil_demand(with_methanol=True)}\n",
    "    | {\"biomass\": get_biomass_demand()}\n",
    "    | {\"process emission\": get_process_emission()},\n",
    "    axis=1,\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2f0fb342-c1b6-4916-ab8b-d0cff9d27475",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_regional_demands(demand_by_region, nodes, \"electricity\", vmax=140)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "38c8c85d-5b22-4b36-8fa8-21625df76024",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_regional_demands(demand_by_region, nodes, \"H2\", vmax=25)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ca7e441a-f183-4b93-867e-0f35aa3a3f04",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_regional_demands(demand_by_region, nodes, \"heat\", vmax=100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d9ed30e5-7a04-469e-b06a-661a66990009",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_regional_demands(\n",
    "    demand_by_region, nodes, \"oil\", series=get_oil_demand(with_methanol=True), vmax=100\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "378d4cd9-431a-4e86-bd25-6e4c67eed732",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_regional_demands(\n",
    "    demand_by_region, nodes, \"solid biomass\", series=get_biomass_demand(), vmax=25\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "50f1ba5f-4230-4f81-a09c-881eb3074437",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_regional_demands(\n",
    "    demand_by_region, nodes, \"process emission\", series=get_process_emission(), vmax=5\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2eae5e22-892c-49d2-87fa-ac427f771e94",
   "metadata": {},
   "outputs": [],
   "source": [
    "if GAS_NETWORK:\n",
    "    plot_regional_demands(demand_by_region, nodes, \"gas\", vmax=16)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b9d950f6-619a-485f-84a9-14ac84b91fa4",
   "metadata": {},
   "source": [
    "## Maps: system-level demands"
   ]
  },
  {
   "cell_type": "raw",
   "id": "00edcc82-c034-4fb4-b6d7-00ea6b1fd67e",
   "metadata": {},
   "source": [
    "def plot_system_demands(df, shape, carrier):\n",
    "    shape[\"value\"] = df.filter(regex=regex[carrier]).sum(axis=1).loc[\"EU\"]\n",
    "\n",
    "    proj = ccrs.EqualEarth()\n",
    "    shape = shape.to_crs(proj.proj4_init)\n",
    "\n",
    "    fig, ax = plt.subplots(figsize=(7, 7), subplot_kw={\"projection\": proj})\n",
    "\n",
    "    shape.plot(\n",
    "        ax=ax,\n",
    "        column=shape.value,\n",
    "        # transform=ccrs.PlateCarree(),\n",
    "        cmap=cmap[carrier],\n",
    "        linewidths=0,\n",
    "        vmin=0,\n",
    "        vmax=shape.value.max() * 2,\n",
    "        legend=True,\n",
    "        legend_kwds={\n",
    "            \"label\": title[carrier],\n",
    "            \"shrink\": 0.7,\n",
    "        },\n",
    "    )\n",
    "\n",
    "    ax.add_feature(cartopy.feature.COASTLINE.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "\n",
    "    ax.set_frame_on(False)\n",
    "    ax.set_facecolor(\"white\")\n",
    "\n",
    "    plt.savefig(f\"{OUTPUT}/demand-map-{carrier}.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "raw",
   "id": "156e53bc-f720-4dae-b543-ca103235e9a9",
   "metadata": {},
   "source": [
    "plot_system_demands(demand_by_region, europe_shape, \"solid biomass\")"
   ]
  },
  {
   "cell_type": "raw",
   "id": "6ed51ba0-096f-4799-a140-12a6784573d3",
   "metadata": {},
   "source": [
    "plot_system_demands(demand_by_region, europe_shape, \"oil\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e20b421a-4601-47ed-b975-7939678d3f34",
   "metadata": {},
   "source": [
    "## Capacity Factors Renewables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1855e884-b19c-4d3b-84a3-88c43c5ee4fb",
   "metadata": {},
   "outputs": [],
   "source": [
    "df = (\n",
    "    n.generators_t.p_max_pu.mean()\n",
    "    .groupby([n.generators.carrier, n.generators.bus.map(n.buses.location)])\n",
    "    .first()\n",
    "    .unstack(0)\n",
    "    .mul(100)\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e136406d-e4b4-4de1-a250-3e2fcdb7bf43",
   "metadata": {},
   "outputs": [],
   "source": [
    "cop = (\n",
    "    n.links_t.efficiency.mean()\n",
    "    .groupby([n.links.carrier, n.links.bus0.map(n.buses.location)])\n",
    "    .first()\n",
    "    .unstack(0)\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "628b9af8-fcb0-4083-bf00-b15220c551f2",
   "metadata": {},
   "outputs": [],
   "source": [
    "report[\"capacity_factor\"] = df.groupby(rename_techs_tyndp, axis=1).mean()\n",
    "report[\"cop\"] = cop.groupby(rename_techs_tyndp, axis=1).mean()\n",
    "\n",
    "df[\"ground-sourced heat pump\"] = cop.filter(like=\"ground\").mean(axis=1)\n",
    "df[\"air-sourced heat pump\"] = cop.filter(like=\"air\").mean(axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "02316ce5-af0e-41dd-9edc-23cf65cb1d34",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_capacity_factors(\n",
    "    df,\n",
    "    geodf,\n",
    "    carrier,\n",
    "    cmap=\"Blues\",\n",
    "    vmax=100,\n",
    "    vmin=0,\n",
    "    label=\"capacity factors [%]\",\n",
    "    fn=None,\n",
    "):\n",
    "    proj = ccrs.EqualEarth()\n",
    "    geodf = geodf.to_crs(proj.proj4_init)\n",
    "\n",
    "    fig, ax = plt.subplots(figsize=(7, 7), subplot_kw={\"projection\": proj})\n",
    "\n",
    "    geodf.plot(\n",
    "        ax=ax,\n",
    "        column=df[carrier].reindex(geodf.index),\n",
    "        cmap=cmap,\n",
    "        linewidths=0,\n",
    "        legend=True,\n",
    "        vmax=vmax,\n",
    "        vmin=vmin,\n",
    "        legend_kwds={\n",
    "            \"label\": label,\n",
    "            \"shrink\": 0.7,\n",
    "            # \"extend\": \"max\",\n",
    "        },\n",
    "    )\n",
    "\n",
    "    ax.add_feature(cartopy.feature.COASTLINE.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "    ax.add_feature(cartopy.feature.BORDERS.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "\n",
    "    ax.set_frame_on(False)\n",
    "    ax.set_facecolor(\"white\")\n",
    "\n",
    "    if fn is None:\n",
    "        plt.savefig(f\"{OUTPUT}/cf-{carrier}.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4de0bcb5-1b2f-4cf8-9451-5428bb430bae",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_capacity_factors(df, nodes, \"onwind\", vmax=55)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7c1c853b-3c9c-42b6-8a31-859d4c748a3f",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_capacity_factors(df, nodes, \"solar\", \"Oranges\", vmax=15)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3ad2395f-0e0a-4bde-88b8-9301b10e52f3",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_capacity_factors(df, offnodes, \"offwind-dc\", vmax=55)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4239d64c-4bea-4fc8-b4b8-760e97822461",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_capacity_factors(df, offnodes, \"offwind-ac\", vmax=55)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8267689c-ebbb-4dc7-bc2c-8052f5b39e61",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_capacity_factors(\n",
    "    df,\n",
    "    nodes,\n",
    "    \"air-sourced heat pump\",\n",
    "    cmap=\"Greens\",\n",
    "    vmax=4,\n",
    "    vmin=2,\n",
    "    label=\"mean COP [-]\",\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "80a6add1-8069-463b-a7c0-359f660611a5",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_capacity_factors(\n",
    "    df,\n",
    "    nodes,\n",
    "    \"ground-sourced heat pump\",\n",
    "    cmap=\"Greens\",\n",
    "    vmax=4,\n",
    "    vmin=2,\n",
    "    label=\"mean COP [-]\",\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ef0c2a30-0df3-41b7-8eb8-b433344389e3",
   "metadata": {},
   "source": [
    "## Demand Totals"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c647dbb2-0c8c-4328-952d-151e308b2197",
   "metadata": {},
   "outputs": [],
   "source": [
    "mapping = {\n",
    "    \"H2 for industry\": \"hydrogen\",\n",
    "    \"H2 for shipping\": \"hydrogen\",\n",
    "    \"agriculture electricity\": \"electricity\",\n",
    "    \"agriculture heat\": \"heat\",\n",
    "    \"agriculture machinery oil\": \"oil\",\n",
    "    \"agriculture machinery oil emissions\": \"emissions\",\n",
    "    \"electricity\": \"electricity\",\n",
    "    \"gas for industry\": \"methane\",\n",
    "    \"industry electricity\": \"electricity\",\n",
    "    \"kerosene for aviation\": \"oil\",\n",
    "    \"land transport EV\": \"electricity\",\n",
    "    \"land transport fuel cell\": \"hydrogen\",\n",
    "    \"low-temperature heat for industry\": \"heat\",\n",
    "    \"naphtha for industry\": \"oil\",\n",
    "    \"shipping methanol\": \"methanol\",\n",
    "    \"shipping methanol emissions\": \"emissions\",\n",
    "    \"oil emissions\": \"emissions\",\n",
    "    \"process emissions\": \"emissions\",\n",
    "    \"residential rural heat\": \"heat\",\n",
    "    \"residential urban decentral heat\": \"heat\",\n",
    "    \"services rural heat\": \"heat\",\n",
    "    \"services urban decentral heat\": \"heat\",\n",
    "    \"solid biomass for industry\": \"solid biomass\",\n",
    "    \"urban central heat\": \"heat\",\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7217bb36-b327-4758-b969-b37ce816cc18",
   "metadata": {},
   "outputs": [],
   "source": [
    "mapping_sector = {\n",
    "    \"H2 for industry\": (\"hydrogen\", \"industry\"),\n",
    "    \"H2 for shipping\": (\"hydrogen\", \"shipping\"),\n",
    "    \"agriculture electricity\": (\"electricity\", \"agriculture\"),\n",
    "    \"agriculture heat\": (\"heat\", \"agriculture\"),\n",
    "    \"agriculture machinery oil\": (\"oil\", \"agriculture\"),\n",
    "    \"agriculture machinery oil emissions\": (\"emissions\", \"agriculture\"),\n",
    "    \"electricity\": (\"electricity\", \"residential\"),\n",
    "    \"gas for industry\": (\"methane\", \"industry\"),\n",
    "    \"industry electricity\": (\"electricity\", \"industry\"),\n",
    "    \"kerosene for aviation\": (\"oil\", \"aviation\"),\n",
    "    \"land transport EV\": (\"electricity\", \"land transport\"),\n",
    "    \"land transport fuel cell\": (\"hydrogen\", \"land transport\"),\n",
    "    \"low-temperature heat for industry\": (\"heat\", \"industry\"),\n",
    "    \"naphtha for industry\": (\"oil\", \"industry\"),\n",
    "    \"shipping methanol\": (\"methanol\", \"shipping\"),\n",
    "    \"shipping methanol emissions\": (\"emissions\", \"other\"),\n",
    "    \"oil emissions\": (\"emissions\", \"other\"),\n",
    "    \"process emissions\": (\"emissions\", \"process\"),\n",
    "    \"residential rural heat\": (\"heat\", \"residential rural\"),\n",
    "    \"residential urban decentral heat\": (\"heat\", \"residential urban\"),\n",
    "    \"services rural heat\": (\"heat\", \"services rural\"),\n",
    "    \"services urban decentral heat\": (\"heat\", \"services urban\"),\n",
    "    \"solid biomass for industry\": (\"solid biomass\", \"industry\"),\n",
    "    \"urban central heat\": (\"heat\", \"district heating\"),\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fe4c65bb-92e8-4936-aa2a-5447036101a2",
   "metadata": {},
   "outputs": [],
   "source": [
    "order = [\n",
    "    \"electricity\",\n",
    "    \"heat\",\n",
    "    \"oil\",\n",
    "    \"methanol\",\n",
    "    \"solid biomass\",\n",
    "    \"methane\",\n",
    "    \"hydrogen\",\n",
    "]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "45824922-5345-4b98-9143-a90c2032b114",
   "metadata": {},
   "outputs": [],
   "source": [
    "df = demand_by_region.sum()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b1e77731-e73f-4457-9ff8-4d7e420ba6da",
   "metadata": {},
   "outputs": [],
   "source": [
    "df.index = pd.MultiIndex.from_tuples([(mapping[i], i) for i in df.index])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b871fe1e-7c33-4a19-b4e4-4b24dd56dc8b",
   "metadata": {},
   "outputs": [],
   "source": [
    "df.drop(\"emissions\", inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8e221e09",
   "metadata": {},
   "outputs": [],
   "source": [
    "colors = config[\"plotting\"][\"tech_colors\"]\n",
    "colors[\"solid biomass\"] = \"seagreen\"\n",
    "colors[\"methane\"] = \"#db6a25\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "56b741c9",
   "metadata": {},
   "outputs": [],
   "source": [
    "pd.DataFrame(df.groupby(level=0).sum().div(1e3).loc[order], columns=[\"\"]).T.plot.bar(\n",
    "    stacked=True,\n",
    "    color=colors,\n",
    "    figsize=(2, 5),\n",
    "    ylim=[0, 20],\n",
    "    ylabel=\"Final energy an non-energy demand [1000 TWh]\",\n",
    ")\n",
    "plt.savefig(f\"{OUTPUT}/total-annual-demand.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "raw",
   "id": "18566d5f-fc8c-404f-b1ff-d5462c8c5668",
   "metadata": {},
   "source": [
    "(n.snapshot_weightings.generators @ loads).div(1e9).sum()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "59300e55-3278-44f7-ba44-e1a5d6cbc14d",
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(5, 5))\n",
    "\n",
    "df.unstack().loc[order].plot.bar(\n",
    "    ax=ax, stacked=True, cmap=\"tab20\", edgecolor=\"k\", ylim=(-100, 4500)\n",
    ")\n",
    "\n",
    "plt.legend(bbox_to_anchor=(1, 1))\n",
    "\n",
    "plt.ylabel(\"Final energy and non-energy demand [TWh/a]\")\n",
    "\n",
    "plt.savefig(f\"{OUTPUT}/demand-by-carrier.pdf\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cb4e8734-2a4c-4c06-b1f9-e9a44cf71f57",
   "metadata": {},
   "outputs": [],
   "source": [
    "df = demand_by_region.sum()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b77f742f-c21b-40ba-b8ec-0d9b71e008d1",
   "metadata": {},
   "outputs": [],
   "source": [
    "df.index = pd.MultiIndex.from_tuples([mapping_sector[i] for i in df.index])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "17ff7dec-e427-4bbe-8399-eb4c40605c20",
   "metadata": {},
   "outputs": [],
   "source": [
    "df = df.loc[df > 0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7ba36311-baa9-4b67-bced-2b642c2811e3",
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(4.5, 4.5))\n",
    "\n",
    "df.unstack().loc[order].T.plot.barh(\n",
    "    ax=ax, color=colors, stacked=True, edgecolor=\"k\", xlim=(-100, 3500)\n",
    ")\n",
    "\n",
    "plt.legend(bbox_to_anchor=(1, 1))\n",
    "\n",
    "plt.xlabel(\"Final energy and non-energy demand [TWh/a]\")\n",
    "\n",
    "plt.savefig(f\"{OUTPUT}/demand-by-sector-carrier.pdf\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4ba1a0d0-3307-4636-b8b5-1dd6d11085f4",
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(4.5, 4.5))\n",
    "\n",
    "df.unstack().loc[order].plot.barh(\n",
    "    ax=ax, stacked=True, cmap=\"tab20\", edgecolor=\"k\", xlim=(-100, 4500)\n",
    ")\n",
    "\n",
    "plt.legend(bbox_to_anchor=(1, 1))\n",
    "\n",
    "plt.xlabel(\"Final energy and non-energy demand [TWh/a]\")\n",
    "\n",
    "plt.savefig(f\"{OUTPUT}/demand-by-carrier-sector.pdf\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1d33b207-8e69-49a8-94cc-c457e186fa23",
   "metadata": {},
   "outputs": [],
   "source": [
    "if not GAS_NETWORK:\n",
    "    plot_system_demands(demand_by_region, europe_shape, \"gas\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "89285f97-86ba-45f5-a6b9-d65e64f7899c",
   "metadata": {
    "tags": []
   },
   "source": [
    "## Demand Time Series"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e8d0d198-19ff-4f5f-a5a2-1bc55dc0d44d",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "kwargs = dict(ylim=[0, 1], xlabel=\"\", figsize=(8, 2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "aceef029-76a7-4e92-b1ac-9cdc80814b73",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "trans = n.loads_t.p_set.filter(like=\"land transport EV\").sum(axis=1)\n",
    "trans /= trans.max()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "875e675b-ed44-4248-84c2-09be04688fc4",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "trans[\"12-2013\"].plot(**kwargs, c=\"#444444\")\n",
    "plt.savefig(f\"{OUTPUT}/bev-demand-month.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "335fe47b-bba9-434f-b61d-7b6f74136b15",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "trans.resample(\"D\").mean().plot(**kwargs, c=\"#444444\")\n",
    "plt.savefig(f\"{OUTPUT}/bev-demand-year.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "636253d4-5b3d-4298-ada1-fbf3cfa2411f",
   "metadata": {},
   "outputs": [],
   "source": [
    "avail = n.links_t.p_max_pu.mean(axis=1)\n",
    "avail[\"12-2013\"].plot(**kwargs, c=\"#baf238\")\n",
    "plt.savefig(f\"{OUTPUT}/bev-availability-month.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1998ac52-cdd6-4a6a-a1ad-98c282fa2444",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "avail.resample(\"D\").mean().plot(**kwargs, c=\"#baf238\")\n",
    "plt.savefig(f\"{OUTPUT}/bev-availability-year.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "64abd051-516c-4397-9137-72c5b8b952a7",
   "metadata": {},
   "outputs": [],
   "source": [
    "heat = n.loads_t.p_set.filter(like=\"heat\").sum(axis=1)\n",
    "heat /= heat.max()\n",
    "heat[\"02-2013\"].plot(**kwargs, c=\"#cc1f1f\")\n",
    "plt.savefig(f\"{OUTPUT}/heat-demand-month.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d15cfd63-797f-4204-acd4-d377fe3df9a0",
   "metadata": {},
   "outputs": [],
   "source": [
    "heat.plot(**kwargs, c=\"#cc1f1f\")\n",
    "plt.savefig(f\"{OUTPUT}/heat-demand-year.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9d508c2c-3d9b-44bd-b90d-06d10978f7e3",
   "metadata": {},
   "outputs": [],
   "source": [
    "elec = n.loads_t.p_set.loc[\n",
    "    :, n.loads.bus.map(n.buses.carrier).isin([\"AC\", \"low voltage\"])\n",
    "].sum(axis=1)\n",
    "elec /= elec.max()\n",
    "elec[\"02-2013\"].plot(**kwargs, c=\"#110d63\")\n",
    "plt.savefig(f\"{OUTPUT}/elec-demand-month.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4642b1b9-59f0-4d8b-9b57-edf506601794",
   "metadata": {},
   "outputs": [],
   "source": [
    "elec.plot(**kwargs, c=\"#110d63\")\n",
    "plt.savefig(f\"{OUTPUT}/elec-demand-year.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "65bf1752",
   "metadata": {},
   "outputs": [],
   "source": [
    "loads = as_dense(n, \"Load\", \"p_set\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "128ede45",
   "metadata": {},
   "outputs": [],
   "source": [
    "heat = loads.filter(like=\"heat\").sum(axis=1).div(1e3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6b3a8de9",
   "metadata": {},
   "outputs": [],
   "source": [
    "elec_nodes = n.loads.bus.map(n.buses.carrier).isin(\n",
    "    [\"AC\", \"low voltage\", \"Li ion\", \"battery\"]\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "207a5d11",
   "metadata": {},
   "outputs": [],
   "source": [
    "elec = loads.loc[:, elec_nodes].sum(axis=1).div(1e3)  # GW"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fca90732",
   "metadata": {},
   "outputs": [],
   "source": [
    "gas_nodes = n.loads.bus.map(n.buses.carrier).isin([\"gas\", \"biogas\", \"gas for industry\"])\n",
    "gas = loads.loc[:, gas_nodes].sum(axis=1).div(1e3)  # GW"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "11673a64",
   "metadata": {},
   "outputs": [],
   "source": [
    "h2_nodes = n.loads.bus.map(n.buses.carrier).isin([\"H2\", \"H2 liquid\"])\n",
    "h2 = loads.loc[:, h2_nodes].sum(axis=1).div(1e3)  # GW"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "328e491d",
   "metadata": {},
   "outputs": [],
   "source": [
    "bio_nodes = n.loads.bus.map(n.buses.carrier).isin(\n",
    "    [\"solid biomass\", \"solid biomass for industry\"]\n",
    ")\n",
    "bio = loads.loc[:, bio_nodes].sum(axis=1).div(1e3)  # GW"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0e871033-fe96-405e-89ee-2e1dae34e0d2",
   "metadata": {},
   "outputs": [],
   "source": [
    "oil_nodes = n.loads.bus.map(n.buses.carrier).isin([\"oil\", \"methanol\"])\n",
    "oil = loads.loc[:, oil_nodes].sum(axis=1).div(1e3)  # GW"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b93e7957",
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(11, 5))\n",
    "\n",
    "gas.plot(ax=ax, color=\"#db6a25\", label=\"Methane\", linestyle=\"--\", linewidth=1)\n",
    "bio.plot(ax=ax, color=\"seagreen\", label=\"Solid Biomass\", linestyle=\"--\", linewidth=1)\n",
    "oil.plot(ax=ax, color=\"k\", label=\"Liquid Hydrocarbons\", linestyle=\"--\", linewidth=1)\n",
    "\n",
    "elec.plot(\n",
    "    ax=ax, color=\"#110d63\", alpha=0.5, label=\"Electricity (3-hourly)\", linewidth=1.5\n",
    ")\n",
    "heat.plot(ax=ax, color=\"#cc1f1f\", alpha=0.5, label=\"Heat (3-hourly)\", linewidth=1.5)\n",
    "h2.plot(ax=ax, color=\"#bf13a0\", alpha=0.5, label=\"Hydrogen (3-hourly)\", linewidth=1.5)\n",
    "\n",
    "\n",
    "elec.rolling(window=\"7D\", center=True).mean().plot(\n",
    "    ax=ax, color=\"#110d63\", label=\"Electricity (7-day rolling mean)\", linewidth=1.5\n",
    ")\n",
    "heat.rolling(window=\"7D\", center=True).mean().plot(\n",
    "    ax=ax, color=\"#cc1f1f\", label=\"Heat (7-day rolling mean)\", linewidth=1.5\n",
    ")\n",
    "h2.rolling(window=\"7D\", center=True).mean().plot(\n",
    "    ax=ax, color=\"#bf13a0\", label=\"Hydrogen (7-day rolling mean)\", linewidth=1.5\n",
    ")\n",
    "\n",
    "plt.xlabel(\"\")\n",
    "plt.legend(ncol=3)\n",
    "plt.ylabel(\"Final energy and non-energy demand [GW]\")\n",
    "\n",
    "plt.savefig(f\"{OUTPUT}/ts-demand.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cb91c5cc-5ffa-456b-bc2e-934981b2567a",
   "metadata": {},
   "source": [
    "## CO2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "76dfb4e2-76a8-479a-a0dc-77476f356243",
   "metadata": {},
   "outputs": [],
   "source": [
    "co2 = pd.read_csv(PATH + \"pypsa-eur-sec/resources/co2_totals.csv\", index_col=0)\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(7, 6))\n",
    "co2.plot.barh(ax=ax, stacked=True, cmap=\"tab20\")\n",
    "plt.savefig(OUTPUT + \"/co2.pdf\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4a469cf6-c3e0-4c9f-b377-ab02ec23f42e",
   "metadata": {},
   "source": [
    "## Industrial Production Fuel Switching"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2dc7ff11-90df-47fd-95cb-0a7f1b8efa8f",
   "metadata": {},
   "outputs": [],
   "source": [
    "fn = PATH + \"pypsa-eur-sec/resources/industrial_production_per_country.csv\"\n",
    "iproduction_today = pd.read_csv(fn, index_col=0).sum()\n",
    "\n",
    "fn = (\n",
    "    PATH + \"pypsa-eur-sec/resources/industrial_production_per_country_tomorrow_2050.csv\"\n",
    ")\n",
    "iproduction_tomorrow = pd.read_csv(fn, index_col=0).sum()\n",
    "\n",
    "fn = PATH + \"pypsa-eur-sec/resources/industry_sector_ratios.csv\"\n",
    "iratios = pd.read_csv(fn, index_col=0)\n",
    "\n",
    "fn = PATH + \"pypsa-eur-sec/resources/industrial_energy_demand_per_country_today.csv\"\n",
    "ienergy_today = (\n",
    "    pd.read_csv(\n",
    "        fn,\n",
    "        index_col=0,\n",
    "        header=[0, 1],\n",
    "    )\n",
    "    .groupby(level=1, axis=1)\n",
    "    .sum()\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "75ec1264-7659-4103-b392-f6efc00ce0bb",
   "metadata": {},
   "outputs": [],
   "source": [
    "ienergy_today.rename({\"gas\": \"methane\"}, inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "05f0c4df-c08a-43f3-b249-1bb945ab43a1",
   "metadata": {},
   "outputs": [],
   "source": [
    "iratios.rename({\"elec\": \"electricity\", \"naphtha\": \"liquid\"}, inplace=True)\n",
    "iratios.loc[\"solid\", :] = iratios.loc[[\"coke\", \"coal\"]].sum()\n",
    "iratios.drop(\n",
    "    [\"coke\", \"coal\", \"process emission\", \"process emission from feedstock\"],\n",
    "    axis=0,\n",
    "    inplace=True,\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "47fc203e-0e72-4f25-9432-9ef401a77540",
   "metadata": {},
   "outputs": [],
   "source": [
    "ienergy_tomorrow = iratios * iproduction_tomorrow / 1e3\n",
    "bc = [\n",
    "    \"Chlorine\",\n",
    "    \"HVC\",\n",
    "    \"HVC (mechanical recycling)\",\n",
    "    \"HVC (chemical recycling)\",\n",
    "    \"Methanol\",\n",
    "]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bd3ccae3-a140-4384-9977-a17899acfdee",
   "metadata": {},
   "outputs": [],
   "source": [
    "ienergy_tomorrow[\"Basic chemicals (without ammonia)\"] = ienergy_tomorrow[bc].sum(axis=1)\n",
    "ienergy_tomorrow.drop(bc, axis=1, inplace=True)\n",
    "ienergy_today[\"DRI + Electric arc\"] = 0.0\n",
    "ienergy_today.loc[\"hydrogen\", :] = 0.0\n",
    "ienergy_tomorrow.loc[\"other\", :] = 0.0\n",
    "ienergy_tomorrow.loc[\"waste\", :] = 0.0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ea02ec5c-f0e9-4124-b60e-90cbc8f75592",
   "metadata": {},
   "outputs": [],
   "source": [
    "ienergy_tomorrow.sort_index(axis=1, inplace=True)\n",
    "ienergy_today.sort_index(axis=1, inplace=True)\n",
    "\n",
    "ienergy_tomorrow.sort_index(axis=0, inplace=True)\n",
    "ienergy_today.sort_index(axis=0, inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a53230c8-4e1b-421e-892f-9e518352ffa3",
   "metadata": {},
   "outputs": [],
   "source": [
    "tech_colors = config[\"plotting\"][\"tech_colors\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6339e15f-6a68-460b-bef1-69ba892c5c6a",
   "metadata": {},
   "outputs": [],
   "source": [
    "tech_colors[\"electricity\"] = \"#ace37f\"\n",
    "tech_colors[\"hydrogen\"] = \"#f073da\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "742fdf9b-5fb1-4ab0-a6ec-81906d09d737",
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(4, 6))\n",
    "\n",
    "ienergy_tomorrow.T.plot.barh(\n",
    "    ax=ax, stacked=True, width=0.3, color=tech_colors, position=1, edgecolor=\"k\"\n",
    ")\n",
    "ienergy_today.T.plot.barh(\n",
    "    ax=ax, stacked=True, width=0.3, color=tech_colors, position=0, edgecolor=\"k\"\n",
    ")\n",
    "\n",
    "\n",
    "plt.xlabel(\"Final energy and non-energy [TWh/a]\")\n",
    "plt.xlim(-50, 1300)\n",
    "plt.ylim(-1, 23)\n",
    "handles, labels = ax.get_legend_handles_labels()\n",
    "n = ienergy_today.shape[0]\n",
    "plt.legend(handles[:n], labels[:n])\n",
    "plt.savefig(OUTPUT + \"/fec_industry_today_tomorrow.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d20118f5-e305-4691-9455-88679b6c418c",
   "metadata": {},
   "outputs": [],
   "source": [
    "df = ienergy_tomorrow.sum(axis=1).sort_values(ascending=False)\n",
    "df = df.loc[df > 0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "01997d5d-f092-4ea8-aae0-802fb33ac240",
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(3, 4))\n",
    "\n",
    "df.plot.bar(ax=ax, width=0.3, edgecolor=\"k\")\n",
    "\n",
    "plt.ylabel(\"Final energy and non-energy [TWh/a]\")\n",
    "plt.xlabel(\"\")\n",
    "\n",
    "plt.savefig(OUTPUT + \"/fec_industry_tomorrow_by_carrier.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e3d23bfd-b585-427c-864f-7419add7b309",
   "metadata": {},
   "source": [
    "## Process Emissions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "95d8441e-d6a0-4ae8-bbbf-dfbf0376b996",
   "metadata": {},
   "outputs": [],
   "source": [
    "pe = (\n",
    "    pd.read_csv(\n",
    "        PATH + \"pypsa-eur-sec/resources/industry_sector_ratios.csv\", index_col=0\n",
    "    )\n",
    "    .filter(like=\"process emission\", axis=0)\n",
    "    .sum()\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f9bc246b-3839-454f-9391-8d54400bdca5",
   "metadata": {},
   "outputs": [],
   "source": [
    "pe_today = iproduction_today * pe / 1e3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "41c3de9f-d735-4e80-bad3-514ae2af008d",
   "metadata": {},
   "outputs": [],
   "source": [
    "pe_tomorrow = iproduction_tomorrow * pe / 1e3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8ed79693-d281-4d04-89d5-ecdfb45ac8c4",
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(4, 6))\n",
    "\n",
    "pe_tomorrow.plot.barh(ax=ax, width=0.3, color=\"peachpuff\", position=1, edgecolor=\"k\")\n",
    "pe_today.plot.barh(ax=ax, width=0.3, color=\"indianred\", position=0, edgecolor=\"k\")\n",
    "\n",
    "plt.xlabel(f\"Process emissions [MtCO$_2$/a]\")\n",
    "plt.xlim(-4, 110)\n",
    "plt.ylim(-1, 27)\n",
    "handles, labels = ax.get_legend_handles_labels()\n",
    "plt.legend(\n",
    "    handles[:n],\n",
    "    [\n",
    "        f\"2050:  {pe_tomorrow.sum():.0f} MtCO$_2$/a\",\n",
    "        f\"today: {pe_today.sum():.0f} MtCO$_2$/a\",\n",
    "    ],\n",
    ")\n",
    "plt.savefig(OUTPUT + \"/process-emissions.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "847cd565-7c7e-4481-98d1-00b56a72861b",
   "metadata": {},
   "source": [
    "## Power Plants"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1e1e80dc-2e68-4c04-bc96-50d3c4816693",
   "metadata": {},
   "outputs": [],
   "source": [
    "df = pd.read_csv(PATH + \"pypsa-eur/resources/powerplants.csv\", index_col=0)\n",
    "\n",
    "df = gpd.GeoDataFrame(df, geometry=gpd.points_from_xy(df.lon, df.lat), crs=\"EPSG:4326\")\n",
    "\n",
    "df.Fueltype = df.Fueltype.str.strip().replace(\"CCGT, Thermal\", \"CCGT\")\n",
    "\n",
    "df = df.loc[\n",
    "    ~df.Fueltype.isin([\"Other\", \"Pv\", \"Storage Technologies\", \"not found\", \"Not Found\"])\n",
    "]\n",
    "\n",
    "colors = {\n",
    "    \"Bioenergy\": \"#80c944\",\n",
    "    \"CCGT\": \"#a85522\",\n",
    "    \"Geothermal\": \"red\",\n",
    "    \"Hard Coal\": \"black\",\n",
    "    \"Lignite\": \"#826837\",\n",
    "    \"Hydro\": \"#235ebc\",\n",
    "    \"Nuclear\": \"#ff8c00\",\n",
    "    \"OCGT\": \"#e0986c\",\n",
    "    \"Oil\": \"#c9c9c9\",\n",
    "    \"Waste\": \"purple\",\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "330d5c80-3848-4222-b046-ed9f59a3227c",
   "metadata": {},
   "outputs": [],
   "source": [
    "crs = ccrs.AlbersEqualArea()\n",
    "\n",
    "df = df.cx[-12:30, 35:72]\n",
    "df = df.to_crs(crs.proj4_init)\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(9.5, 9.5), subplot_kw={\"projection\": crs})\n",
    "\n",
    "ax.add_feature(cartopy.feature.COASTLINE.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "ax.add_feature(cartopy.feature.BORDERS.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "\n",
    "df.plot(\n",
    "    ax=ax,\n",
    "    column=\"Fueltype\",\n",
    "    markersize=df[\"Capacity\"] / 35,\n",
    "    alpha=0.75,\n",
    "    legend=True,\n",
    "    cmap=mcolors.ListedColormap(\n",
    "        pd.Series(df.Fueltype.unique()).sort_values().map(colors).values\n",
    "    ),\n",
    "    legend_kwds=dict(title=\"Technology (size ~ capacity)\", frameon=False),\n",
    ")\n",
    "\n",
    "\n",
    "ax.set_frame_on(False)\n",
    "ax.set_facecolor(\"white\")\n",
    "plt.savefig(OUTPUT + \"/powerplants.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ea648e6f-e467-43a1-8051-f8a0608d03e0",
   "metadata": {},
   "source": [
    "## Hotmaps Industrial Sites"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3562badc-1db4-4a40-84bd-13f489c13c71",
   "metadata": {},
   "outputs": [],
   "source": [
    "def prepare_hotmaps_database():\n",
    "    \"\"\"\n",
    "    Load hotmaps database of industrial sites.\n",
    "    \"\"\"\n",
    "\n",
    "    df = pd.read_csv(\n",
    "        PATH + \"pypsa-eur-sec/data/Industrial_Database.csv\", sep=\";\", index_col=0\n",
    "    )\n",
    "\n",
    "    df[[\"srid\", \"coordinates\"]] = df.geom.str.split(\";\", expand=True)\n",
    "\n",
    "    # remove those sites without valid locations\n",
    "    df.drop(df.index[df.coordinates.isna()], inplace=True)\n",
    "\n",
    "    df[\"coordinates\"] = gpd.GeoSeries.from_wkt(df[\"coordinates\"])\n",
    "\n",
    "    gdf = gpd.GeoDataFrame(df, geometry=\"coordinates\", crs=\"EPSG:4326\")\n",
    "\n",
    "    return gdf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f5c0c7f0-9aca-4179-864e-ed64430f042c",
   "metadata": {},
   "outputs": [],
   "source": [
    "crs = ccrs.AlbersEqualArea()\n",
    "\n",
    "hotmaps = prepare_hotmaps_database()\n",
    "hotmaps = hotmaps.cx[-12:30, 35:72]\n",
    "hotmaps = hotmaps.to_crs(crs.proj4_init)\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(9.5, 9.5), subplot_kw={\"projection\": crs})\n",
    "\n",
    "ax.add_feature(cartopy.feature.COASTLINE.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "ax.add_feature(cartopy.feature.BORDERS.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "\n",
    "hotmaps.plot(\n",
    "    ax=ax,\n",
    "    column=\"Subsector\",\n",
    "    markersize=hotmaps[\"Emissions_ETS_2014\"] / 3e4,\n",
    "    alpha=0.5,\n",
    "    legend=True,\n",
    "    # transform=ccrs.epsg(3395),\n",
    "    legend_kwds=dict(title=\"Industry Sector (size ~ emissions)\", frameon=False),\n",
    ")\n",
    "\n",
    "\n",
    "ax.set_frame_on(False)\n",
    "ax.set_facecolor(\"white\")\n",
    "\n",
    "plt.savefig(OUTPUT + \"/hotmaps.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ae392a84-c879-4297-baa4-18e48e140670",
   "metadata": {},
   "source": [
    "## Salt Caverns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ef057c15-2b1f-43b1-9e47-2f757ebb898c",
   "metadata": {},
   "outputs": [],
   "source": [
    "cavern_nodes = pd.read_csv(\n",
    "    f\"{PATH}/pypsa-eur-sec/resources/salt_cavern_potentials_s_{CLUSTERS}.csv\",\n",
    "    index_col=0,\n",
    ")\n",
    "cavern_nodes = cavern_nodes.where(cavern_nodes > 0.5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "938db9f0-c10e-454d-9b98-ae825ef49d69",
   "metadata": {},
   "outputs": [],
   "source": [
    "report[\"salt_caverns\"] = cavern_nodes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "34aea126-6a37-4d73-9a33-1c5f56f8a09d",
   "metadata": {},
   "outputs": [],
   "source": [
    "cavern_regions = gpd.read_file(\n",
    "    f\"{PATH}pypsa-eur/resources/regions_onshore_elec_s_{CLUSTERS}.geojson\"\n",
    ").set_index(\"name\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "75380497-1217-41a8-9bb1-a26b43410b0e",
   "metadata": {},
   "outputs": [],
   "source": [
    "cavern_offregions = gpd.read_file(\n",
    "    f\"{PATH}/pypsa-eur/resources/regions_offshore_elec_s_{CLUSTERS}.geojson\"\n",
    ").set_index(\"name\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2bfbb98e-41d0-43fc-81c0-e6c3bf2b4e12",
   "metadata": {},
   "outputs": [],
   "source": [
    "cavern_regions.to_crs(crs.proj4_init).total_bounds"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "619cd68d-1749-45fb-b979-c039df2e480e",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_salt_caverns_by_node(\n",
    "    cavern_nodes,\n",
    "    cavern_regions,\n",
    "    storage_type=\"onshore\",\n",
    "    cmap=\"GnBu\",\n",
    "    vmin=1,\n",
    "    vmax=3000,\n",
    "    fn=None,\n",
    "    label=r\"H$_2$ Storage Potential [TWh]\",\n",
    "):\n",
    "    crs = ccrs.EqualEarth()\n",
    "\n",
    "    cavern_regions = cavern_regions.to_crs(crs.proj4_init)\n",
    "\n",
    "    fig, ax = plt.subplots(figsize=(7, 7), subplot_kw={\"projection\": crs})\n",
    "\n",
    "    cavern_regions.plot(\n",
    "        ax=ax,\n",
    "        column=cavern_nodes[storage_type].reindex(cavern_regions.index),\n",
    "        # transform=ccrs.PlateCarree(),\n",
    "        cmap=cmap,\n",
    "        linewidths=0,\n",
    "        legend=True,\n",
    "        vmin=vmin,\n",
    "        vmax=vmax,\n",
    "        legend_kwds={\n",
    "            \"label\": label,\n",
    "            \"shrink\": 0.7,\n",
    "            \"extend\": \"max\",\n",
    "        },\n",
    "        norm=mcolors.LogNorm(vmin=1, vmax=vmax),\n",
    "    )\n",
    "\n",
    "    plt.title(f\"{storage_type.capitalize()} Salt Cavern H$_2$ Storage Potentials\")\n",
    "\n",
    "    plt.xlim(-1.2e6, 2.6e6)\n",
    "    plt.ylim(4.3e6, 7.8e6)\n",
    "\n",
    "    # plt.xlim(-14.5, 31.5)\n",
    "    # plt.ylim(34, 72)\n",
    "\n",
    "    ax.add_feature(cartopy.feature.COASTLINE.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "    ax.add_feature(cartopy.feature.BORDERS.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "\n",
    "    ax.set_frame_on(False)\n",
    "    ax.set_facecolor(\"white\")\n",
    "\n",
    "    if fn is None:\n",
    "        plt.savefig(\n",
    "            f\"{OUTPUT}/cavern-potentials-{storage_type}.pdf\", bbox_inches=\"tight\"\n",
    "        )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ac7049b1-6c16-4c81-9a46-2bd42734b14c",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_salt_caverns_by_node(cavern_nodes, cavern_regions, storage_type=\"onshore\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4bc4f0d1-4d26-4f64-ae25-487b0e5b5e57",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_salt_caverns_by_node(cavern_nodes, cavern_regions, storage_type=\"nearshore\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ba1a0b95-7fd8-43e0-a93c-ff915b8d1366",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_salt_caverns_by_node(cavern_nodes, cavern_offregions, storage_type=\"offshore\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d2376e96-d3b7-478d-b1ec-bd798315bc99",
   "metadata": {},
   "outputs": [],
   "source": [
    "caverns = gpd.read_file(\n",
    "    f\"{PATH}/pypsa-eur-sec/data/h2_salt_caverns_GWh_per_sqkm.geojson\"\n",
    ")\n",
    "\n",
    "crs = ccrs.EqualEarth()\n",
    "\n",
    "caverns = caverns.to_crs(crs.proj4_init)\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(7, 7), subplot_kw={\"projection\": crs})\n",
    "\n",
    "ax.add_feature(cartopy.feature.COASTLINE.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "ax.add_feature(cartopy.feature.BORDERS.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "\n",
    "caverns.plot(\n",
    "    ax=ax,\n",
    "    column=\"storage_type\",\n",
    "    cmap=\"tab10_r\",  # \"tab10_r\",\n",
    "    legend=True,\n",
    "    linewidth=0,\n",
    "    # transform=ccrs.epsg(3395),\n",
    "    legend_kwds=dict(\n",
    "        title=\"Salt Caverns for\\nHydrogen Storage\", frameon=False, loc=(0.2, 0.8)\n",
    "    ),\n",
    ")\n",
    "\n",
    "\n",
    "plt.xlim(-1e6, 2.6e6)\n",
    "plt.ylim(4.3e6, 7.8e6)\n",
    "\n",
    "ax.set_frame_on(False)\n",
    "ax.set_facecolor(\"white\")\n",
    "\n",
    "plt.savefig(OUTPUT + \"/caverns.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1751b0e3-b587-4864-a04b-3c6d549e37de",
   "metadata": {},
   "source": [
    "## Biomass potentials"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a975fc61-8393-4d36-aee5-db4cc923ea70",
   "metadata": {},
   "outputs": [],
   "source": [
    "fn = PATH + f\"pypsa-eur-sec/resources/biomass_potentials_s_{CLUSTERS}.csv\"\n",
    "bio = pd.read_csv(fn, index_col=0).div(1e6).drop(\"not included\", axis=1)  # TWh/a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ab8a8289-8a68-4918-9e8b-986627ecfa2b",
   "metadata": {},
   "outputs": [],
   "source": [
    "report[\"biomass_potentials\"] = bio"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "44abd552-ae0a-43a9-97e2-ea908fd29541",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_biomass_potentials(bio, nodes, kind, fn=None):\n",
    "    crs = ccrs.EqualEarth()\n",
    "    nodes = nodes.to_crs(crs.proj4_init)\n",
    "\n",
    "    fig, ax = plt.subplots(figsize=(7, 7), subplot_kw={\"projection\": crs})\n",
    "\n",
    "    ax.add_feature(cartopy.feature.COASTLINE.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "    ax.add_feature(cartopy.feature.BORDERS.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "\n",
    "    nkind = \"disregarded biomass\" if kind == \"not included\" else kind\n",
    "    label = f\"{nkind} potentials [TWh/a]\"\n",
    "\n",
    "    nodes.plot(\n",
    "        ax=ax,\n",
    "        column=bio[kind],\n",
    "        cmap=\"Greens\",  # \"tab10_r\",\n",
    "        legend=True,\n",
    "        linewidth=0,\n",
    "        legend_kwds={\n",
    "            \"label\": label,\n",
    "            \"shrink\": 0.7,\n",
    "            \"extend\": \"max\",\n",
    "        },\n",
    "    )\n",
    "\n",
    "    total = bio[kind].sum()\n",
    "\n",
    "    ax.text(-0.8e6, 7.4e6, f\"total: {total:.0f} TWh/a\", fontsize=15, color=\"#343434\")\n",
    "\n",
    "    plt.xlim(-1e6, 2.6e6)\n",
    "    plt.ylim(4.3e6, 7.8e6)\n",
    "\n",
    "    ax.set_frame_on(False)\n",
    "    ax.set_facecolor(\"white\")\n",
    "\n",
    "    if fn is None:\n",
    "        plt.savefig(OUTPUT + f\"/biomass-{kind}.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3da2f01f-6a65-438c-a491-2471b807d474",
   "metadata": {},
   "outputs": [],
   "source": [
    "for kind in bio.columns:\n",
    "    plot_biomass_potentials(bio, nodes, kind)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6ea2de59-1819-4696-ae17-2d517ca7ed01",
   "metadata": {},
   "source": [
    "## Ammonia Production"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bb918a6b-657b-4b59-8dca-f4864bff68e0",
   "metadata": {},
   "outputs": [],
   "source": [
    "fn = PATH + f\"pypsa-eur-sec/resources/ammonia_production.csv\"\n",
    "df = pd.read_csv(fn, index_col=0)  # kt/a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fb8000b6-c5bc-4be3-97dd-94d24cff719a",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_ammonia_production(df, cts, year, fn=None):\n",
    "    year = str(year)\n",
    "\n",
    "    crs = ccrs.EqualEarth()\n",
    "    cts = cts.to_crs(crs.proj4_init)\n",
    "\n",
    "    fig, ax = plt.subplots(figsize=(7, 7), subplot_kw={\"projection\": crs})\n",
    "\n",
    "    ax.add_feature(cartopy.feature.COASTLINE.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "    ax.add_feature(cartopy.feature.BORDERS.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "\n",
    "    cts.plot(\n",
    "        ax=ax,\n",
    "        column=df[year].reindex(cts.index),\n",
    "        cmap=\"Purples\",  # \"tab10_r\",\n",
    "        legend=True,\n",
    "        linewidth=0,\n",
    "        legend_kwds={\n",
    "            \"label\": \"Ammonia Production [kt]\",\n",
    "            \"shrink\": 0.7,\n",
    "        },\n",
    "    )\n",
    "\n",
    "    total = df[year].sum()\n",
    "    ax.text(\n",
    "        -0.85e6, 7.4e6, f\"total: {total:.0f} kt ({year})\", fontsize=12, color=\"#343434\"\n",
    "    )\n",
    "\n",
    "    plt.xlim(-1e6, 2.6e6)\n",
    "    plt.ylim(4.3e6, 7.8e6)\n",
    "\n",
    "    ax.set_frame_on(False)\n",
    "\n",
    "    if fn is None:\n",
    "        plt.savefig(OUTPUT + f\"/ammonia-{year}.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4c79177d-ca13-4e49-a3ef-b78d1d04fb45",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_ammonia_production(df, cts, 2017)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "910f2473-c216-4aef-9a6b-36f9bb0d8dc2",
   "metadata": {},
   "source": [
    "## Gas Network Unclustered"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "02a30a3b-9ee0-4011-931b-b43fe0b303d5",
   "metadata": {},
   "outputs": [],
   "source": [
    "lng_fn = PATH + \"pypsa-eur-sec/data/gas_network/scigrid-gas/data/IGGIELGN_LNGs.geojson\"\n",
    "entry_fn = (\n",
    "    PATH\n",
    "    + \"pypsa-eur-sec/data/gas_network/scigrid-gas/data/IGGIELGN_BorderPoints.geojson\"\n",
    ")\n",
    "prod_fn = (\n",
    "    PATH\n",
    "    + \"pypsa-eur-sec/data/gas_network/scigrid-gas/data/IGGIELGN_Productions.geojson\"\n",
    ")\n",
    "planned_lng_fn = PATH + \"pypsa-eur-sec/data/gas_network/planned_LNGs.csv\"\n",
    "\n",
    "regions_onshore = PATH + \"pypsa-eur/resources/regions_onshore_elec_s_181.geojson\"\n",
    "regions_offshore = PATH + \"pypsa-eur/resources/regions_offshore_elec_s_181.geojson\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ce516e26-e2be-40aa-8b1f-5f0bb3415926",
   "metadata": {},
   "outputs": [],
   "source": [
    "regions = load_bus_regions(regions_onshore, regions_offshore)\n",
    "\n",
    "# add a buffer to eastern countries because some\n",
    "# entry points are still in Russian or Ukrainian territory.\n",
    "buffer = 9000  # meters\n",
    "eastern_countries = [\"FI\", \"EE\", \"LT\", \"LV\", \"PL\", \"SK\", \"HU\", \"RO\"]\n",
    "add_buffer_b = regions.index.str[:2].isin(eastern_countries)\n",
    "regions.loc[add_buffer_b] = (\n",
    "    regions[add_buffer_b].to_crs(3035).buffer(buffer).to_crs(4326)\n",
    ")\n",
    "\n",
    "countries = regions.index.str[:2].unique().str.replace(\"GB\", \"UK\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "050835d8-a4eb-402d-8692-069c67de8129",
   "metadata": {},
   "outputs": [],
   "source": [
    "pts = build_gas_input_locations(\n",
    "    lng_fn,\n",
    "    planned_lng_fn,\n",
    "    entry_fn,\n",
    "    prod_fn,\n",
    "    countries,\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f881862a-383e-4296-9d1d-04ed762bf0ad",
   "metadata": {},
   "outputs": [],
   "source": [
    "sums = pts.groupby(\"type\").p_nom.sum()\n",
    "\n",
    "for t in [\"lng\", \"production\", \"pipeline\"]:\n",
    "    pts.loc[pts[\"type\"] == t, \"p_nom\"] /= sums[t]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "290ae1c4-8879-45aa-b4f8-229869b05f9d",
   "metadata": {},
   "outputs": [],
   "source": [
    "pts[\"type\"] = pts[\"type\"].replace(\n",
    "    dict(production=\"Fossil Extraction\", lng=\"LNG Terminal\", pipeline=\"Entrypoint\")\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d4738e3b-a45e-4494-8f12-aa2a3951f4c6",
   "metadata": {},
   "outputs": [],
   "source": [
    "fn = PATH + \"pypsa-eur-sec/resources/gas_network.csv\"\n",
    "df = pd.read_csv(fn, index_col=0)\n",
    "for col in [\"geometry\"]:\n",
    "    df[col] = df[col].apply(wkt.loads)\n",
    "\n",
    "df = gpd.GeoDataFrame(df, geometry=\"geometry\", crs=\"EPSG:4326\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "92ddd252-a764-4cda-9248-52a6db862708",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "crs = ccrs.AlbersEqualArea()\n",
    "\n",
    "df = df.to_crs(crs.proj4_init)\n",
    "pts = pts.to_crs(crs.proj4_init)\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(9.5, 9.5), subplot_kw={\"projection\": crs})\n",
    "\n",
    "ax.add_feature(cartopy.feature.COASTLINE.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "ax.add_feature(cartopy.feature.BORDERS.with_scale(\"50m\"), linewidth=0.2, zorder=2)\n",
    "\n",
    "df.plot(\n",
    "    ax=ax,\n",
    "    column=df[\"p_nom\"].div(1e3),\n",
    "    linewidths=df[\"p_nom\"].clip(upper=50e3).div(2e4),\n",
    "    cmap=\"Spectral_r\",\n",
    "    vmax=60,\n",
    "    legend=True,\n",
    "    legend_kwds=dict(label=\"Gas Pipeline Capacity [GW]\", shrink=0.5, extend=\"max\"),\n",
    ")\n",
    "\n",
    "pts.plot(\n",
    "    ax=ax,\n",
    "    column=\"type\",\n",
    "    markersize=pts[\"p_nom\"] * 300,\n",
    "    legend=True,\n",
    ")\n",
    "\n",
    "\n",
    "ax.set_frame_on(False)\n",
    "ax.set_facecolor(\"white\")\n",
    "plt.savefig(OUTPUT + \"/gas_network.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7797c025-43d3-4bc0-9689-347caf5e56c1",
   "metadata": {},
   "source": [
    "## Report CSV"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2d98b4dd-755b-422e-a20b-14d9b4f7edf6",
   "metadata": {},
   "outputs": [],
   "source": [
    "report = {k: v.round(2) for k, v in report.items()}\n",
    "report = pd.concat(report, axis=1)\n",
    "\n",
    "report.to_csv(f\"{OUTPUT}/report.csv\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a4d98807-7562-4aa8-b4fe-9f7afbbbe18d",
   "metadata": {},
   "outputs": [],
   "source": [
    "ds = report.stack(level=1)\n",
    "ds.index.names = [\"bus\", \"carrier\"]\n",
    "ds.to_xarray().to_netcdf(\"resources.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "991e4527-76ad-42bf-bb4c-ba08499e5f76",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "494d35a45045e3fb9015f1d2cd094f8e1540ef3ea8b346f6f51a6ceb6be203c2"
  },
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}