scripts/OSM/1a_initialize_osm.ipynb

{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "source": [
    "<img src=\"../../images/BikeDNA_logo.svg\" width=\"250\"  alt=\"BikeDNA logo\" style=\"display:block; margin-left: auto; margin-right: auto;\">\n",
    "<a href=\"https://github.com/anerv/BikeDNA\">Github</a>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 1a. Initialize OSM data\n",
    "\n",
    "This notebook:\n",
    "\n",
    "- Loads the polygon defining the study area and then creates a grid overlay for the study area.\n",
    "- Downloads street network data for the study area using OSMnx.\n",
    "- Creates a network only with bicycle infrastructure (with queries defined in `config.yml`).\n",
    "- Creates additional attributes in the data to be used in the analysis."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "source": [
    "<div class=\"alert alert-block alert-info\">\n",
    "<b>Prerequisites &amp; Input/Output</b>\n",
    "<p>\n",
    "<span style=\"font-family:courier;\">config.yml</span> must be set in advance.\n",
    "    \n",
    "Output files of this notebook are saved to the <span style=\"font-family:courier;\">../data/OSM/[study_area]/processed/</span> subfolder and used in notebooks 1b, 3a and 3b.\n",
    "\n",
    "When the network is created from the OSM data, each edge in the simplified and non-simplified network is assigned a unique edge ID (<i>not</i> OSM ID). To relate the edge ID back to the OSM ID, the dataframes with network edges, <span style=\"font-family:courier;\">osm_edges</span> and <span style=\"font-family:courier;\">osm_edges_simplified</span>, can be used.\n",
    "</p>\n",
    "</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Sections**\n",
    "* [Load data for study area and create analysis grid](#Load-data-for-study-area-and-create-analysis-grid)\n",
    "* [Download and preprocess OSM data](#Download-and-preprocess-OSM-data)\n",
    "\n",
    "<br />"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load data for study area and create analysis grid\n",
    "\n",
    "This step:\n",
    "\n",
    "- Loads settings for the analysis from the configuration file `config.yml`.\n",
    "- Reads data for the study area.\n",
    "- Creates a grid overlay of the study area, with grid cell size as defined in `config.yml`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "source": [
    "<div class=\"alert alert-block alert-warning\">\n",
    "<b>Troubleshooting</b>\n",
    "<p>\n",
    "If in the cell below the functions located in the <i>bikedna</i> folder do not load properly: Check that <span style=\"font-family:courier;\">pip install -e .</span> was run successfully (see instructions in <span style=\"font-family:courier;\">README.md</span>).\n",
    "</p>\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "outputs": [],
   "source": [
    "# Load libraries, settings and data\n",
    "\n",
    "import json\n",
    "import pickle\n",
    "import warnings\n",
    "from datetime import datetime\n",
    "\n",
    "import contextily as cx\n",
    "import geopandas as gpd\n",
    "import matplotlib.pyplot as plt\n",
    "import networkx as nx\n",
    "import osmnx as ox\n",
    "import numpy as np\n",
    "import yaml\n",
    "\n",
    "from bikedna import evaluation_functions as eval_func\n",
    "from bikedna import graph_functions as graph_func\n",
    "from bikedna import simplification_functions as simp_func\n",
    "from bikedna import plotting_functions as plot_func\n",
    "\n",
    "%run ../settings/plotting.py\n",
    "%run ../settings/yaml_variables.py\n",
    "%run ../settings/paths.py\n",
    "\n",
    "warnings.filterwarnings(\"ignore\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Load data for study area\n",
    "\n",
    "The study area is defined by the user-provided polygon. It will be used for the computation of **global** results, i.e. quality metrics based on all data in the study area."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Load polygon for study area\n",
    "\n",
    "study_area_poly = gpd.read_file(study_area_poly_fp)\n",
    "\n",
    "if study_area_poly.crs == None:\n",
    "    print(\"Please assign a crs to the study area polygon!\")\n",
    "\n",
    "if study_area_poly.crs != study_crs:\n",
    "    study_area_poly = study_area_poly.to_crs(study_crs)\n",
    "\n",
    "assert study_area_poly.crs == study_crs\n",
    "\n",
    "area = study_area_poly.area.values[0]\n",
    "print(f\"The size of the study area is {area / 1000000:.2f} km2.\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "outputs": [],
   "source": [
    "# Plot study area\n",
    "\n",
    "set_renderer(renderer_map)\n",
    "print(\"This study area will be used for the analysis:\")\n",
    "\n",
    "fig, ax = plt.subplots(1, figsize=pdict[\"fsmap\"])\n",
    "study_area_poly.plot(ax=ax, edgecolor=pdict[\"base\"], facecolor='none',alpha=pdict[\"alpha_back\"],hatch='/')\n",
    "ax.set_axis_off()\n",
    "cx.add_basemap(ax, crs=study_area_poly.crs, source=cx_tile_1)\n",
    "ax.set_title(f\"{area_name} study area ({area / 1000000:.2f} sqkm)\")\n",
    "\n",
    "plot_func.save_fig(fig, osm_results_static_maps_fp + \"area_size_osm\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Create analysis grid"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "source": [
    "<div class=\"alert alert-block alert-info\">\n",
    "<b>User configurations</b>\n",
    "<p>\n",
    "The grid is used to compute <b>local</b> results for many of the quality metrics, i.e. on grid cell level. The cell size of the grid cells must be specified in the configuration file. Smaller grid cells will increase the computation time but give a more granular result. The grid cell size should not be smaller than the scale at which you would expect actual differences in data quality.\n",
    "</p>\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create and plot grid\n",
    "\n",
    "set_renderer(renderer_map)\n",
    "\n",
    "grid = eval_func.create_grid_geometry(study_area_poly, grid_cell_size)\n",
    "grid[\"grid_id\"] = grid.index\n",
    "\n",
    "fig, ax = plt.subplots(1, figsize=pdict[\"fsmap\"])\n",
    "grid_cell_area = (grid_cell_size / 1000) ** 2  # in km2\n",
    "print(f\"The grid contains {len(grid)} square cells with a side length of {grid_cell_size} m and an area of {grid_cell_area:.2f} km2.\")\n",
    "print(\"This grid will be used for local (grid cell level) analysis:\")\n",
    "grid.plot(ax=ax, facecolor=\"none\", edgecolor=pdict[\"base\"], alpha=pdict[\"alpha_back\"],linewidth=0.5)\n",
    "ax.set_axis_off()\n",
    "cx.add_basemap(ax, crs=study_area_poly.crs, source=cx_tile_1)\n",
    "ax.set_title(\n",
    "    f\"{area_name} study area ({len(grid)} grid cells, side length {grid_cell_size}m)\"\n",
    ")\n",
    "\n",
    "plot_func.save_fig(fig, osm_results_static_maps_fp + \"area_grid_osm\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Download and preprocess OSM data\n",
    "\n",
    "This step:\n",
    "\n",
    "- Downloads data from OpenStreetMap using OSMnx.\n",
    "- Projects the data to the chosen CRS.\n",
    "- Creates a subnetwork consisting only of bicycle infrastructure.\n",
    "- Classifies all edges in the bicycle network based on whether they are protected or unprotected bicycle infrastructure, how they have been digitized, and whether they allow for bidirectional travel or not.\n",
    "- Simplifies the network.\n",
    "- Creates copies of all edge and node data sets indexed by their intersecting grid cell."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**OSM data model**\n",
    "\n",
    "In OSM, street network data are stored using *nodes* (points) and *ways* (lines). In BikeDNA, OSM data are converted to a network structure consisting of *nodes* and *edges* (we use the terminology used in OSMnx). Edges represents the actual infrastructure, such as bike lanes and paths, while nodes represents the start and end points for the edges, as well as all intersections. For further details, read more about the [OSM data model](https://alga.win.tue.nl/tutorials/openstreetmap/) and the [network data model](https://networkx.org/documentation/stable/reference/introduction.html)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "source": [
    "<details><summary><b>Network simplification</b></summary>\n",
    "\n",
    "When converting a set of linestrings to network (graph) format, not all vertices (nodes) are of equal meaning. For the geometry of edges, all nodes are used to define the shape of the edge. For the topology of the network, however, only those nodes that are endpoints or intersection points with other edges are needed, while all other (so-called 'interstitial') nodes do not add any information. To compare the structure and true ratio between nodes and edges in a network, a simplified network representation which only includes nodes at endpoints and intersections, or where the value of important attributes changes, is required. An additional advantage of simplifying the network is the resulting substantial reduction of the number of nodes and edges, which makes computational routines much faster.\n",
    "\n",
    "Note that the simplification function only simplifies and merges edges that have identical values for important attributes such as geometry type, bicycle direction, and protection level. \n",
    "\n",
    "<p align=\"center\">\n",
    "<img src='../../images/network_simplification_illustration.png' width=300/>\n",
    "\n",
    "*Non-simplified network (left) and simplified network (right)*.\n",
    "\n",
    "</p>\n",
    "\n",
    "As part of the simplification routine, in cases where there are several edges between the same pair of nodes ('parallel edges' or 'multiedges'), only one of the edges is retained.\n",
    "\n",
    "To read more about the modified OSMnx simplification [(Boeing, 2017)](https://geoffboeing.com/publications/osmnx-complex-street-networks/) used here, we refer to this [GitHub repository](https://github.com/anerv/NERDS_osmnx) which contains both the simplification functions, explanation of the logic and a demonstration.\n",
    "\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "source": [
    "<div class=\"alert alert-block alert-danger\">\n",
    "<b>Warning</b>\n",
    "<p>\n",
    "Depending on the size of your study area, getting data through OSMnx can take a little while (for a city around the size of Copenhagen it can take two minutes).\n",
    "</p>\n",
    "</div>\n",
    "\n",
    "<div class=\"alert alert-block alert-info\">\n",
    "<b>User configurations</b>\n",
    "<p>\n",
    "In the following step, the dedicated bicycle infrastructure is extracted from the general street network. Depending on the type of bicycle infrastructure and OSM tagging tradition in your area, the queries used for retrieving the bicycle network might need to be adapted (this is done in the configuration file under <i>bicycle_infrastructure_queries</i>).\n",
    "</p>\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "outputs": [],
   "source": [
    "# Download OSM data\n",
    "\n",
    "ox.settings.use_cache = False\n",
    "ox.settings.useful_tags_way = osm_way_tags\n",
    "\n",
    "graph_osm = ox.graph_from_polygon(\n",
    "    study_area_poly.to_crs(\"EPSG:4326\").loc[0, \"geometry\"],\n",
    "    network_type=\"bike\",\n",
    "    simplify=False,\n",
    "    retain_all=True,\n",
    "    truncate_by_edge=False,\n",
    "    clean_periphery=True,\n",
    ")\n",
    "\n",
    "print(\"OSM data retrieved successfully!\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "outputs": [],
   "source": [
    "# Convert into graph\n",
    "\n",
    "# Project graph to chosen crs\n",
    "graph_osm = ox.project_graph(graph_osm, to_crs=study_crs)\n",
    "\n",
    "# Get osm_edges and osm_nodes\n",
    "osm_nodes, osm_edges = ox.graph_to_gdfs(graph_osm)\n",
    "\n",
    "# Clean up column names (removing ':' and upper case letters)\n",
    "osm_edges = graph_func.clean_col_names(osm_edges)\n",
    "osm_nodes = graph_func.clean_col_names(osm_nodes)\n",
    "\n",
    "# Add attribute on whether bicycle infra exist or not (to be used by e.g. simplification function)\n",
    "osm_edges[\"bicycle_infrastructure\"] = \"no\"\n",
    "\n",
    "for q in bicycle_infrastructure_queries.values():\n",
    "\n",
    "    try:\n",
    "        ox_filtered = osm_edges.query(q)\n",
    "\n",
    "    except Exception:\n",
    "        print(\"Exception occured when quering with:\", q)\n",
    "        print(\"Please check if the columns used in the query are present in the data\")\n",
    "\n",
    "    osm_edges.loc[ox_filtered.index, \"bicycle_infrastructure\"] = \"yes\"\n",
    "\n",
    "osm_edges.bicycle_infrastructure.value_counts()\n",
    "\n",
    "bicycle_edges = osm_edges.loc[osm_edges.bicycle_infrastructure == \"yes\"].copy()\n",
    "\n",
    "bicycle_graph = graph_osm.edge_subgraph(bicycle_edges.index)\n",
    "print(\"Successfully created network only with bicycle infrastructure!\")\n",
    "print(\"\\n\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Augment and classify edge types\n",
    "\n",
    "assert len(bicycle_edges) == len(bicycle_graph.edges)\n",
    "\n",
    "bicycle_graph = graph_func.update_key_values(bicycle_graph)\n",
    "\n",
    "bicycle_nodes, bicycle_edges = ox.graph_to_gdfs(bicycle_graph)\n",
    "bicycle_edges = graph_func.clean_col_names(bicycle_edges)\n",
    "bicycle_nodes = graph_func.clean_col_names(bicycle_nodes)\n",
    "\n",
    "bicycle_edges[\"bicycle_infrastructure\"] = \"yes\"\n",
    "bicycle_infrastructure_dict = bicycle_edges[\"bicycle_infrastructure\"].to_dict()\n",
    "nx.set_edge_attributes(\n",
    "    bicycle_graph, bicycle_infrastructure_dict, \"bicycle_infrastructure\"\n",
    ")\n",
    "\n",
    "# Add \"multiple edge\" attribute to each edge of the graph\n",
    "for e in bicycle_graph.edges:\n",
    "    bicycle_graph.edges[e][\"multiedge\"] = 0\n",
    "\n",
    "# Count multiple edges and store info at key 0\n",
    "# (simplification routine keeps only key 0)\n",
    "for u, v, k in bicycle_graph.edges:\n",
    "    if k > 0:\n",
    "        bicycle_graph.edges[u, v, 0][\"multiedge\"] += 1\n",
    "\n",
    "# Define whether bicycle infrastructure is in both or one direction and whether it is one- or bidirectional\n",
    "bicycle_edges = eval_func.simplify_bicycle_tags(bicycle_edges,simplify_tags_queries)\n",
    "\n",
    "for key, value in bicycle_edges.bicycle_bidirectional.value_counts().items():\n",
    "    perc = np.round(100*value/len(bicycle_edges),2)\n",
    "    print(f\"Edges where 'bicycle_bidirectional' is {key}: {value} out of {len(bicycle_edges)} ({perc}%)\")\n",
    "print(\"\\n\")\n",
    "\n",
    "for key, value in bicycle_edges.bicycle_geometries.value_counts().items():\n",
    "    perc = np.round(100*value/len(bicycle_edges),2)\n",
    "    print(f\"Edges where the geometry type is '{key}': {value} out of {len(bicycle_edges)} ({perc}%)\")\n",
    "print(\"\\n\")\n",
    "\n",
    "bicycle_bidirectional_dict = bicycle_edges[\"bicycle_bidirectional\"].to_dict()\n",
    "nx.set_edge_attributes(\n",
    "    bicycle_graph, bicycle_bidirectional_dict, \"bicycle_bidirectional\"\n",
    ")\n",
    "\n",
    "bicycle_geometries_dict = bicycle_edges[\"bicycle_geometries\"].to_dict()\n",
    "nx.set_edge_attributes(bicycle_graph, bicycle_geometries_dict, \"bicycle_geometries\")\n",
    "\n",
    "# Classify edges as protected or unprotected\n",
    "bicycle_edges = eval_func.define_protected_unprotected(\n",
    "    bicycle_edges, osm_bicycle_infrastructure_type\n",
    ")\n",
    "\n",
    "# Set edges attributes for column protected\n",
    "bicycle_protected_dict = bicycle_edges[\"protected\"].to_dict()\n",
    "nx.set_edge_attributes(bicycle_graph, bicycle_protected_dict, \"protected\")\n",
    "\n",
    "for key, value in bicycle_edges.protected.value_counts().items():\n",
    "    perc = np.round(100*value/len(bicycle_edges),2)\n",
    "    print(f\"Edges where the protection level is '{key}': {value} out of {len(bicycle_edges)} ({perc}%)\")\n",
    "print(\"\\n\")\n",
    "\n",
    "# Simplify bicycle network\n",
    "bicycle_graph_simplified = simp_func.simplify_graph(\n",
    "    bicycle_graph,\n",
    "    attributes=[\n",
    "        \"bicycle_infrastructure\",\n",
    "        \"bicycle_bidirectional\",\n",
    "        \"bicycle_geometries\",\n",
    "        \"protected\",\n",
    "    ],\n",
    "    remove_rings=False,\n",
    ")\n",
    "\n",
    "# Get undirected network representation\n",
    "bicycle_graph_simplified = ox.get_undirected(bicycle_graph_simplified)\n",
    "bicycle_graph = ox.get_undirected(bicycle_graph)\n",
    "\n",
    "bicycle_nodes_simplified, bicycle_edges_simplified = ox.graph_to_gdfs(\n",
    "    bicycle_graph_simplified, nodes=True\n",
    ")\n",
    "\n",
    "# Measure the length of bicycle infrastructure (!=from length of edges)\n",
    "bicycle_edges_simplified[\"length\"] = bicycle_edges_simplified.geometry.length\n",
    "bicycle_edges_simplified[\"infrastructure_length\"] = bicycle_edges_simplified.apply(\n",
    "    lambda x: eval_func.measure_infrastructure_length(\n",
    "        edge=x.geometry,\n",
    "        geometry_type=x.bicycle_geometries,\n",
    "        bidirectional=x.bicycle_bidirectional,\n",
    "        bicycle_infrastructure=x.bicycle_infrastructure,\n",
    "    ),\n",
    "    axis=1,\n",
    ")\n",
    "\n",
    "# Creating columns with unique feature id (required by some functions)\n",
    "bicycle_edges[\"edge_id\"] = bicycle_edges.reset_index().index\n",
    "bicycle_edges_simplified[\"edge_id\"] = bicycle_edges_simplified.reset_index().index\n",
    "\n",
    "osm_id_dict = bicycle_edges[\"edge_id\"].to_dict()\n",
    "nx.set_edge_attributes(bicycle_graph, osm_id_dict, \"edge_id\")\n",
    "\n",
    "osm_id_dict_s = bicycle_edges_simplified[\"edge_id\"].to_dict()\n",
    "nx.set_edge_attributes(bicycle_graph_simplified, osm_id_dict_s, \"edge_id\")\n",
    "\n",
    "# Saving data to graphs\n",
    "osm_l_dict = osm_edges[\"length\"].to_dict()\n",
    "nx.set_edge_attributes(bicycle_graph, osm_l_dict, \"length\")\n",
    "\n",
    "osm_l_dict_s = bicycle_edges_simplified[\"length\"].to_dict()\n",
    "nx.set_edge_attributes(bicycle_graph_simplified, osm_l_dict_s, \"length\")\n",
    "\n",
    "osm_il_dict = bicycle_edges_simplified[\"infrastructure_length\"].to_dict()\n",
    "nx.set_edge_attributes(bicycle_graph_simplified, osm_il_dict, \"infrastructure_length\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# Plot global statistics (by network length)\n",
    "\n",
    "set_renderer(renderer_plot)\n",
    "ymax = bicycle_edges_simplified.infrastructure_length.sum() * 1.2 /1000 # conversion m to km\n",
    "subplotdata = []\n",
    "subplotbar_labels = []\n",
    "subplottitle = []\n",
    "subplotx_positions = []\n",
    "\n",
    "# attribute: bidirectional\n",
    "x = []\n",
    "y = []\n",
    "\n",
    "for value in np.unique(bicycle_edges.bicycle_bidirectional):\n",
    "    x.append(str(value))\n",
    "    y.append(bicycle_edges[bicycle_edges[\"bicycle_bidirectional\"]==value].length.sum()/1000) # conversion m to km\n",
    "\n",
    "plot_func.make_bar_plot(\n",
    "    data = y,\n",
    "    bar_labels = [\"unidirectional\", \"bidirectional\"],\n",
    "    y_label = \"Network edge length [km]\",\n",
    "    x_positions = [i for i in range(len(x))],\n",
    "    title = f\"{area_name}\\n OSM edges by bidirectionality\",\n",
    "    bar_colors = len(x)*[pdict[\"osm_base\"]],\n",
    "    filepath = osm_results_plots_fp + \"edge_global_bidirectional\",\n",
    "    ylim = ymax,\n",
    "    figsize=pdict[\"fsbar_small\"]\n",
    "); plt.close()\n",
    "subplotdata.append(tuple(y))\n",
    "subplotbar_labels.append([\"unidirectional\", \"bidirectional\"])\n",
    "subplottitle.append(f\"\\n OSM edges by bidirectionality\")\n",
    "subplotx_positions.append([i for i in range(len(x))])\n",
    "\n",
    "# attribute: bicycle geometry (mapping type)\n",
    "x = []\n",
    "y = []\n",
    "\n",
    "for value in np.unique(bicycle_edges.bicycle_geometries):\n",
    "    x.append(value)\n",
    "    y.append(bicycle_edges[bicycle_edges[\"bicycle_geometries\"]==value].length.sum()/1000) # conversion m to km\n",
    "\n",
    "plot_func.make_bar_plot(\n",
    "    data = y,\n",
    "    bar_labels = x,\n",
    "    y_label = \"Network edge lengths [km]\",\n",
    "    x_positions = [i for i in range(len(x))],\n",
    "    title = f\"{area_name}\\n OSM edges by mapping type\",\n",
    "    bar_colors = len(x)*[pdict[\"osm_base\"]],\n",
    "    filepath = osm_results_plots_fp + \"edge_global_mappingtype\",\n",
    "    ylim = ymax,\n",
    "    figsize=pdict[\"fsbar_small\"]\n",
    "); plt.close()\n",
    "subplotdata.append(tuple(y))\n",
    "subplotbar_labels.append(x)\n",
    "subplottitle.append(f\"{area_name}\\n OSM edges by mapping type\")\n",
    "subplotx_positions.append([i for i in range(len(x))])\n",
    "\n",
    "# attribute: protection\n",
    "x = []\n",
    "y = []\n",
    "\n",
    "for value in np.unique(bicycle_edges.protected):\n",
    "    x.append(value)\n",
    "    y.append(bicycle_edges[bicycle_edges[\"protected\"]==value].length.sum()/1000) # conversion m to km\n",
    "\n",
    "plot_func.make_bar_plot(\n",
    "    data = y,\n",
    "    bar_labels = x,\n",
    "    y_label = \"Network edge lengths [km]\",\n",
    "    x_positions = [i for i in range(len(x))],\n",
    "    title = f\"{area_name}\\n OSM edges by infrastructure type\",\n",
    "    bar_colors = len(x)*[pdict[\"osm_base\"]],\n",
    "    filepath = osm_results_plots_fp + \"edge_global_infrastructuretype\",\n",
    "    ylim = ymax,\n",
    "    figsize=pdict[\"fsbar_small\"]\n",
    "); plt.close()\n",
    "subplotdata.append(tuple(y))\n",
    "subplotbar_labels.append(x)\n",
    "subplottitle.append(f\"\\n OSM edges by infrastructure type\")\n",
    "subplotx_positions.append([i for i in range(len(x))])\n",
    "\n",
    "# Next to each other\n",
    "plot_func.make_bar_subplots(\n",
    "    subplot_data=subplotdata,\n",
    "    nrows=1,\n",
    "    ncols=3,\n",
    "    bar_labels=subplotbar_labels,\n",
    "    y_label=[\"Network edge lengths [km]\", \"\", \"\"],\n",
    "    x_positions=subplotx_positions,\n",
    "    title=subplottitle,\n",
    "    bar_colors=len(x)*[pdict[\"osm_base\"]],\n",
    "    filepath=osm_results_plots_fp + \"edge_global\",\n",
    "    ylim=ymax\n",
    ");\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Overview of data from OSM\n",
    "graph_area = eval_func.get_graph_area(bicycle_nodes, study_area_poly, study_crs)\n",
    "print(f\"The network covers an area of {graph_area/ 1000000:.2f} km2.\")\n",
    "print(\n",
    "    f\"The length of the OSM network with bicycle infrastructure is {bicycle_edges_simplified.infrastructure_length.sum()/1000 :.2f} km.\"\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot network\n",
    "\n",
    "set_renderer(renderer_map)\n",
    "fig, ax = plt.subplots(1, figsize=pdict[\"fsmap\"])\n",
    "bicycle_edges_simplified.plot(ax=ax, color=pdict[\"osm_base\"], linewidth=pdict[\"line_base\"])\n",
    "bicycle_nodes_simplified.plot(ax=ax, color=pdict[\"osm_emp\"], markersize=pdict[\"mark_emp\"])\n",
    "study_area_poly.plot(ax=ax, edgecolor=pdict[\"base\"], facecolor=\"None\", linewidth=1)\n",
    "cx.add_basemap(ax=ax, crs=study_area_poly.crs, source=cx_tile_2)\n",
    "ax.set_axis_off()\n",
    "ax.set_title(f\"{area_name}, OSM network\")\n",
    "\n",
    "plot_func.save_fig(fig, osm_results_static_maps_fp + \"area_network_osm\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "outputs": [],
   "source": [
    "# Joining grid cell id to network elements\n",
    "\n",
    "osm_edges_simp_joined = gpd.overlay(\n",
    "    bicycle_edges_simplified.reset_index(), grid, how=\"intersection\"\n",
    ")\n",
    "osm_edges_joined = gpd.overlay(bicycle_edges.reset_index(), grid, how=\"intersection\")\n",
    "\n",
    "osm_nodes_simp_joined = gpd.overlay(\n",
    "    bicycle_nodes_simplified.reset_index(), grid, how=\"intersection\"\n",
    ")\n",
    "osm_nodes_joined = gpd.overlay(bicycle_nodes.reset_index(), grid, how=\"intersection\")\n",
    "\n",
    "# Count features in each grid cell\n",
    "all_data_osm = [\n",
    "    osm_edges_joined,\n",
    "    osm_nodes_joined,\n",
    "    osm_edges_simp_joined,\n",
    "    osm_nodes_simp_joined,\n",
    "]\n",
    "labels_osm = [\"osm_edges\", \"osm_nodes\", \"osm_simplified_edges\", \"osm_simplified_nodes\"]\n",
    "\n",
    "for data, label in zip(all_data_osm, labels_osm):\n",
    "\n",
    "    df = eval_func.count_features_in_grid(data, label)\n",
    "\n",
    "    grid = eval_func.merge_results(grid, df, \"left\")\n",
    "\n",
    "# Recomputing infrastructure length for data joined to grid\n",
    "osm_edges_simp_joined[\"length\"] = osm_edges_simp_joined[\"geometry\"].length\n",
    "osm_edges_simp_joined[\"infrastructure_length\"] = osm_edges_simp_joined.apply(\n",
    "    lambda x: eval_func.measure_infrastructure_length(\n",
    "        edge=x.geometry,\n",
    "        geometry_type=x.bicycle_geometries,\n",
    "        bidirectional=x.bicycle_bidirectional,\n",
    "        bicycle_infrastructure=x.bicycle_infrastructure,\n",
    "    ),\n",
    "    axis=1,\n",
    ")\n",
    "\n",
    "assert round(bicycle_edges_simplified.infrastructure_length.sum() / 10000, 0) == round(\n",
    "    osm_edges_simp_joined.infrastructure_length.sum() / 10000, 0\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "outputs": [],
   "source": [
    "# Save OSM data\n",
    "exec(open(\"../settings/save_osmdata.py\").read())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "source": [
    "***"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "outputs": [],
   "source": [
    "from time import strftime\n",
    "print(\"Time of analysis: \" + strftime(\"%a, %d %b %Y %H:%M:%S\"))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "celltoolbar": "Tags",
  "kernelspec": {
   "display_name": "bikedna",
   "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.11.0"
  },
  "properties": {
   "title": {
    "description": "Load and process OSM data",
    "type": "string"
   }
  },
  "vscode": {
   "interpreter": {
    "hash": "25d89d047ff73a4f08e1f58df0505313a978ba3f16552d4fb4e98a48d36b765b"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}