BRNO_Data_process_analysis.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import geopandas as gpd\n",
    "import pandas as pd\n",
    "import contextily as ctx\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from shapely.geometry import Polygon\n",
    "import seaborn as sns\n",
    "from mpl_toolkits.axes_grid1 import make_axes_locatable\n",
    "from scipy.interpolate import griddata\n",
    "from matplotlib import colors"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# PART 1: DATA PREPROCESSING"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Dataset: Surface Temperature"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf_st_clean = gpd.read_file('Surface_temperature-shp.zip')\n",
    "gdf_st_clean"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(12, 12))\n",
    "ax.set_axis_off()\n",
    "\n",
    "divider = make_axes_locatable(ax)\n",
    "cax = divider.append_axes(\"right\", size=\"2%\", pad=0.2)\n",
    "gdf_st_clean.to_crs(epsg=3857).plot(ax=ax, markersize=10, alpha=0.6, column='Teplota', cmap='plasma', legend=True, cax=cax)\n",
    "cb_ax = fig.axes[1]\n",
    "cb_ax.tick_params(labelsize=12)\n",
    "cb_ax.set_ylabel(\"\\nTeplota [°C]\", fontsize=15)\n",
    "\n",
    "ctx.add_basemap(ax, crs=gdf_st_clean.to_crs(epsg=3857).crs.to_string(), source=ctx.providers.Stamen.Terrain)\n",
    "\n",
    "fig.savefig('01_Surface_Temperature_clean.png', bbox_inches='tight', dpi=250)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# GeoDataFrame Surface Temperature\n",
    "gdf_st = gdf_st_clean.copy()\n",
    "gdf_st = gdf_st.to_crs(epsg=3857)\n",
    "gdf_st = gdf_st.set_geometry(gdf_st.centroid)\n",
    "gdf_st"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Dataset: Building height"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf_bh_clean = gpd.read_file('Building_heights-shp.zip')\n",
    "gdf_bh_clean"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(12, 12))\n",
    "ax.set_axis_off()\n",
    "\n",
    "gdf_bh_clean.to_crs(epsg=3857).plot(ax=ax)\n",
    "ctx.add_basemap(ax, crs=gdf_st_clean.to_crs(epsg=3857).crs.to_string(), source=ctx.providers.Stamen.Terrain)\n",
    "\n",
    "fig.savefig('02_Building_Height_clean.png', bbox_inches='tight', dpi=250)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf_bh = gdf_bh_clean.copy().to_crs(epsg=5514)  # Krovak projection to calculate area\n",
    "gdf_bh = gdf_bh[gdf_bh['vyska_komp'].notna()]   # drop NaN values\n",
    "gdf_bh['area'] = gdf_bh.area                    # calculate area in m2\n",
    "gdf_bh = gdf_bh.to_crs(epsg=3857)               # back to Mercator projection\n",
    "gdf_bh"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "x_min, y_min, x_max, y_max = gdf_bh.total_bounds\n",
    "\n",
    "# Split the bh map into NxN grid\n",
    "N = 100\n",
    "x_side = (x_max - x_min)/N\n",
    "y_side = (y_max - y_min)/N\n",
    "\n",
    "grid = []\n",
    "for j in range(N):\n",
    "    for i in range(N):\n",
    "        x_point_list = [x_min + i*x_side, x_min + (i+1)*x_side, x_min + (i+1)*x_side, x_min + i*x_side, x_min + i*x_side]\n",
    "        y_point_list = [y_min + j*y_side, y_min + j*y_side, y_min + (j+1)*y_side, y_min + (j+1)*y_side, y_min + j*y_side]\n",
    "        grid.append(Polygon(zip(x_point_list, y_point_list)))\n",
    "\n",
    "gdf_grid = gpd.GeoDataFrame(grid, columns=['geometry'])\n",
    "gdf_grid = gdf_grid.set_crs('epsg:3857')\n",
    "gdf_grid"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Sides of 1 tile in grid\n",
    "print(\"x dimension of 1 tile: %.2f m\" %x_side)\n",
    "print(\"y dimension of 1 tile: %.2f m\" %y_side)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(12, 12))\n",
    "ax.set_axis_off()\n",
    "gdf_grid[\"geometry\"].boundary.plot(color='k', ax=ax, linewidth=0.5, alpha=0.4)\n",
    "gdf_bh.to_crs(epsg=3857).plot(ax=ax)\n",
    "ctx.add_basemap(ax, crs=gdf_bh.crs.to_string(), source=ctx.providers.Stamen.Terrain)\n",
    "\n",
    "fig.savefig('03_Building_Height_grid.png', bbox_inches='tight', dpi=250)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Create/Load dataframe of intersections between grid and building height\n",
    "To create dataframe of intersections uncomment commented lines below. It takes approx 180 min.\n",
    "Otherwise dataframe will be loaded from pickled file."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# gdf_bh_grid = gpd.GeoDataFrame(columns=['geometry','rectangle_ID', 'building_ID'])\n",
    "\n",
    "# for ID, rectangle in enumerate(gdf_grid.geometry):                          # for each tile in grid\n",
    "#     df_tmp = gpd.GeoDataFrame(columns=['geometry','rectangle_ID'])          # prepare empty GeoDataFrame\n",
    "#     df_tmp.geometry = gdf_bh.intersection(rectangle)                        # get intersection of all buildings with current tile\n",
    "#     df_tmp = df_tmp[~df_tmp.is_empty]                                       # drop all empty intersections\n",
    "#     df_tmp[\"rectangle_ID\"] = ID                                             # mark tile ID\n",
    "#     df_tmp = df_tmp.reset_index().rename(columns={'index': 'building_ID'})\n",
    "#     all_df = all_df.append(df_tmp)\n",
    "#     if ID%1000==0:\n",
    "#         print(ID)\n",
    "# gdf_bh_grid.to_pickle(\"./gdf_bh_grid.pkl\")\n",
    "\n",
    "# Load dataframe already created\n",
    "gdf_bh_grid = pd.read_pickle(\"./gdf_bh_grid.pkl\")\n",
    "gdf_bh_grid = gdf_bh_grid.set_crs('epsg:3857')\n",
    "gdf_bh_grid['area_building_in'] = gdf_bh_grid.to_crs(epsg=5514).area    # Krovak to calculate area of parts of buildings inside squares\n",
    "gdf_bh_grid = gdf_bh_grid.to_crs(epsg=3857)                             # Back to Mercator for nicer visualization\n",
    "gdf_bh_grid"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf_bh_grid = gdf_bh_grid.merge(gdf_bh, left_on='building_ID', right_index=True)\n",
    "gdf_bh_grid = gdf_bh_grid.assign(typ_budovy=gdf_bh_grid[\"zpusob_vyu\"].values) # make copy of values from \"zpusob_vyu\" column into \"building_type\" because they will be used twice\n",
    "gdf_bh_grid"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_bh_grid = gdf_bh_grid.pivot_table(index='rectangle_ID', columns=['typ_budovy'],\n",
    "                                values=['zpusob_vyu','vyska_komp','area_building_in', \"OBJECTID_1\"],\n",
    "                                aggfunc={\"zpusob_vyu\":\"count\",\"vyska_komp\":np.mean, \"area_building_in\":np.sum, \"OBJECTID_1\":np.mean}\n",
    "                                ).rename(columns={'zpusob_vyu': 'cnt', 'vyska_komp': 'avg_height', 'area_building_in': 'b_area_in_rect', \"OBJECTID_1\":\"b_area_ratio\"})\n",
    "df_bh_grid['b_area_ratio'] = np.nan\n",
    "df_bh_grid['b_area_ratio'] = df_bh_grid['b_area_in_rect'].div(df_bh_grid['b_area_in_rect'].sum(axis=1), axis=0)\n",
    "df_bh_grid"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Height of tiles\n",
    "Calculate Weighted Mean of buildings height according to the area they occupy in given tile"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf_grid_extended = gdf_grid.loc[df_bh_grid.index.tolist(), :].copy()\n",
    "gdf_grid_extended[\"WM_Height\"] = np.nan\n",
    "for idx in gdf_grid_extended.index:\n",
    "    area = df_bh_grid.loc[idx, 'b_area_in_rect'].loc[~df_bh_grid.loc[idx, 'b_area_in_rect'].isnull()].values\n",
    "    height = df_bh_grid.loc[idx, 'avg_height'].loc[~df_bh_grid.loc[idx, 'avg_height'].isnull()].values\n",
    "    weighted_mean = np.dot(area,height) / np.sum(area)\n",
    "    gdf_grid_extended.loc[idx, 'WM_Height'] = weighted_mean\n",
    "gdf_grid_extended"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(12, 15))\n",
    "ax.set_axis_off()\n",
    "\n",
    "divider = make_axes_locatable(ax)\n",
    "cax = divider.append_axes(\"right\", size=\"2%\", pad=0.2)\n",
    "gdf_grid_extended.plot(ax=ax, alpha=0.8, column='WM_Height', cmap='viridis', legend=True, cax=cax)\n",
    "cb_ax = fig.axes[1]\n",
    "cb_ax.tick_params(labelsize=12)\n",
    "cb_ax.set_ylabel(\"\\nVýška [m]\", fontsize=15)\n",
    "\n",
    "ctx.add_basemap(ax, crs=gdf_st_clean.to_crs(epsg=3857).crs.to_string(), source=ctx.providers.Stamen.Terrain)\n",
    "plt.tight_layout()\n",
    "\n",
    "fig.savefig('04_Building_Height.png', bbox_inches='tight', dpi=250)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Temperature of tiles\n",
    "Take temperature points that are inside building height dataset with 10% border, to not interpolate the whole temperature map when buildings are only in some areas"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Outer boundary of grid\n",
    "gdf_1Square = gpd.GeoDataFrame()\n",
    "\n",
    "x_point_list = [x_min, x_min + (N+1)*x_side, x_min + (N+1)*x_side, x_min, x_min]\n",
    "y_point_list = [y_min, y_min, y_min + (N+1)*y_side, y_min + (N+1)*y_side, y_min]\n",
    "\n",
    "gdf_1Square = gpd.GeoDataFrame(geometry=[Polygon(zip(x_point_list, y_point_list))])\n",
    "gdf_1Square = gdf_1Square.set_crs('epsg:3857')\n",
    "gdf_1Square"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(12, 12))\n",
    "ax.set_axis_off()\n",
    "gdf_bh.plot(ax=ax)\n",
    "gdf_1Square[\"geometry\"].boundary.plot(color='k', ax=ax, linewidth=0.5)\n",
    "gdf_1Square.buffer(0.1*(x_max-x_min)).boundary.plot(color='k', ax=ax, linewidth=0.5)\n",
    "ctx.add_basemap(ax, crs=gdf_1Square.crs.to_string(), alpha=0.8)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Get temperature of points inside buffer\n",
    "gdf_TiB = gpd.GeoDataFrame()\n",
    "gdf_TiB.geometry = gdf_st.intersection(gdf_1Square.buffer(0.1*(x_max-x_min)).geometry.iloc[0])\n",
    "gdf_TiB = gdf_TiB.set_crs(\"epsg:3857\")\n",
    "\n",
    "gdf_TiB = gdf_TiB[~gdf_TiB.geometry.is_empty]                                   # Remove empty points\n",
    "gdf_TiB = gdf_TiB.merge(gdf_st[\"Teplota\"], left_index=True, right_index=True)   # Add temperature values\n",
    "gdf_TiB"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(12, 15))\n",
    "ax.set_axis_off()\n",
    "divider = make_axes_locatable(ax)\n",
    "cax = divider.append_axes(\"right\", size=\"2%\", pad=0.2)\n",
    "\n",
    "gdf_grid[\"geometry\"].boundary.plot(color='k', ax=ax, linewidth=0.5, alpha=0.2)\n",
    "gdf_1Square.buffer(0.1*(x_max-x_min)).boundary.plot(color='k', ax=ax, alpha=0.2, linewidth=0.5)\n",
    "gdf_TiB.plot(ax=ax, markersize=10, alpha=0.7, column='Teplota', cmap='plasma', legend=True, cax=cax)\n",
    "cb_ax = fig.axes[1]\n",
    "cb_ax.tick_params(labelsize=12)\n",
    "cb_ax.set_ylabel(\"\\nTeplota [°C]\", fontsize=15)\n",
    "\n",
    "ctx.add_basemap(ax, crs=gdf_st_clean.to_crs(epsg=3857).crs.to_string(), source=ctx.providers.Stamen.Terrain)\n",
    "plt.tight_layout()\n",
    "\n",
    "fig.savefig('05_Surface_Temperature_centroids.png', bbox_inches='tight', dpi=250)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf_grid[\"Temperature\"] = np.nan\n",
    "\n",
    "points = np.vstack((gdf_TiB.centroid.x, gdf_TiB.centroid.y)).T\n",
    "values = gdf_TiB.Teplota.values\n",
    "\n",
    "grid_x, grid_y = np.meshgrid(gdf_grid.loc[df_bh_grid.index.tolist(),:].centroid.x.values, gdf_grid.loc[df_bh_grid.index.tolist(),:].centroid.y.values)\n",
    "grid_z = griddata(points, values, (grid_x, grid_y), method='linear')\n",
    "gdf_grid.loc[df_bh_grid.index.tolist(), \"Temperature\"] = np.diagonal(grid_z)\n",
    "gdf_grid"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(12, 15))\n",
    "ax.set_axis_off()\n",
    "divider = make_axes_locatable(ax)\n",
    "cax = divider.append_axes(\"right\", size=\"2%\", pad=0.2)\n",
    "\n",
    "gdf_grid.plot(ax=ax, markersize=10, alpha=0.7, column='Temperature', cmap='plasma', legend=True, cax=cax)\n",
    "cb_ax = fig.axes[1]\n",
    "cb_ax.tick_params(labelsize=12)\n",
    "cb_ax.set_ylabel(\"\\nTeplota [°C]\", fontsize=15)\n",
    "\n",
    "ctx.add_basemap(ax, crs=gdf_st_clean.to_crs(epsg=3857).crs.to_string(), source=ctx.providers.Stamen.Terrain)\n",
    "fig.savefig('06_Surafe_Temperature_grid.png', bbox_inches='tight', dpi=250)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Add temperature to the df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_bh_grid.loc[:, \"Temperature\"] = gdf_grid.loc[df_bh_grid.index.tolist(), \"Temperature\"]\n",
    "df_bh_grid"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# PART 2: ANALYSIS"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A \"quick look\" at Building Height dataset to better understand and see the count of each building type, average height, average area and total area covered"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf_bh['area2'] = gdf_bh['area']\n",
    "gdf_bh.groupby('zpusob_vyu').agg({\"OBJECTID_1\":\"count\",\"vyska_komp\":np.mean, \"area\":np.mean, \"area2\":np.sum}\n",
    "                                ).rename(columns={'OBJECTID_1': 'cnt', 'vyska_komp': 'avg_height', 'area': 'avg_area', 'area2': 'sum_area'}\n",
    "                                ).sort_values(by='cnt', ascending=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Split tiles into 3 groups based on temperature\n",
    "### Tiles with:\n",
    "- Temperature below 21°C\n",
    "- Temperature between 21°C and 25°C\n",
    "- Temperature above 25°C"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "labels = ['nad 25°C', '21°C ÷ 25°C', 'pod 21°C']\n",
    "values = [df_bh_grid[df_bh_grid[\"Temperature\"] > 25 ]['b_area_ratio'].sum(axis=0).sum(axis=0) * 100 / df_bh_grid['b_area_ratio'].sum(axis=0).sum(axis=0),\n",
    "        df_bh_grid[(df_bh_grid[\"Temperature\"] <= 25) & (df_bh_grid[\"Temperature\"] >= 21)]['b_area_ratio'].sum(axis=0).sum(axis=0) * 100 / df_bh_grid['b_area_ratio'].sum(axis=0).sum(axis=0),\n",
    "        df_bh_grid[df_bh_grid[\"Temperature\"] < 21 ]['b_area_ratio'].sum(axis=0).sum(axis=0) * 100 / df_bh_grid['b_area_ratio'].sum(axis=0).sum(axis=0)]\n",
    "explode = (0.15,0.05,0.1)\n",
    "colors = [sns.color_palette('plasma')[i] for i in [5,2,0]]\n",
    "\n",
    "fig1, ax1 = plt.subplots(figsize=(5,3))\n",
    "patches, labels, pct_texts = ax1.pie(values, labels=labels, labeldistance=1.05, autopct='%1.1f%%', pctdistance=0.7,\n",
    "                                explode=explode, colors=colors, shadow=False, startangle=39, rotatelabels =True, textprops={'fontsize': 12})\n",
    "\n",
    "for label, pct_text in zip(labels, pct_texts):\n",
    "    pct_text.set_rotation(label.get_rotation())\n",
    "\n",
    "for i, (label, pct_text) in enumerate(zip(labels, pct_texts)):\n",
    "    label.set_rotation(0)\n",
    "    if i in [1,2]:\n",
    "        pct_text.set_color('white')\n",
    "\n",
    "ax1.axis('equal')\n",
    "\n",
    "plt.axis(\"equal\")\n",
    "plt.tight_layout()\n",
    "fig1.savefig('07_Pie_temperature_intervals.svg', pad_inches=0, transparent=True, bbox_inches='tight', dpi=250)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Percentage of the tiles (with at least 1 building) inside temperature intervals<br>\n",
    "- 26.2% tiles T < 21°C <br>\n",
    "- 68.8% tiles 25°C >= T >= 21°C <br>\n",
    "-  5.0% tiles T > 25°C"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(8, 8))\n",
    "ax.set_axis_off()\n",
    "\n",
    "gdf_grid[gdf_grid['Temperature']>25].plot(ax=ax, markersize=10, alpha=0.8,  color=sns.color_palette('plasma')[5])\n",
    "gdf_grid[(gdf_grid['Temperature']<=25) & (gdf_grid['Temperature'] >= 21)].plot(ax=ax, markersize=10, alpha=0.8, color=sns.color_palette('plasma')[2])\n",
    "gdf_grid[gdf_grid['Temperature']<21].plot(ax=ax, markersize=10, alpha=0.8, color=sns.color_palette('plasma')[0])\n",
    "\n",
    "ctx.add_basemap(ax, crs=gdf_st_clean.to_crs(epsg=3857).crs.to_string(), source=ctx.providers.Stamen.Terrain)\n",
    "plt.tight_layout()\n",
    "fig.savefig('08_Surafe_Temperature_grid_intervals.png', bbox_inches='tight', dpi=250)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Buildings in tiles with temperature below 21°C"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_below21 = pd.Series(df_bh_grid[df_bh_grid[\"Temperature\"] < 21 ]['b_area_ratio'].sum(axis=0) * 100 / df_bh_grid[df_bh_grid[\"Temperature\"] < 21 ]['b_area_ratio'].sum(axis=0).sum(axis=0), name=\"Occupy\").to_frame()\n",
    "df_below21 = df_below21[df_below21.values>2].sort_values(by='Occupy', ascending=False)\n",
    "df_below21"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note: sum of area percentage in \"Occupy\" column is not 100% because buildings with area < 2% have been dropped"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Average height calculated from weighted mean height values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "avg_height_below21 = gdf_grid_extended.loc[df_bh_grid[df_bh_grid[\"Temperature\"] < 21 ].index, \"WM_Height\"].sum() / len(df_bh_grid[df_bh_grid[\"Temperature\"] < 21 ].index)\n",
    "print(\"Average height of tiles with temperature below 21°C is %.2f m.\" %avg_height_below21)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "labels = df_below21.index.values\n",
    "values = df_below21.values.flatten()\n",
    "explode = np.size(df_below21.values)*(0.02,)\n",
    "explode = list(explode)\n",
    "explode[0] = 0.05\n",
    "explode = tuple(explode)\n",
    "\n",
    "colors = plt.cm.plasma(np.linspace(0.12,0.28,np.size(df_below21.values)))\n",
    "fig1, ax1 = plt.subplots(figsize=(8,5))\n",
    "patches, labels, pct_texts = ax1.pie(values, labels=labels, labeldistance=1.05, autopct='%1.1f%%', pctdistance=0.8,\n",
    "                                explode=explode, colors=colors, shadow=False, startangle=39, rotatelabels =True, textprops={'fontsize': 12})\n",
    "\n",
    "for label, pct_text in zip(labels, pct_texts):\n",
    "    pct_text.set_rotation(label.get_rotation())\n",
    "    pct_text.set_color('white')\n",
    "\n",
    "for label, pct_text in zip(labels, pct_texts):\n",
    "    label.set_rotation(0)\n",
    "\n",
    "ax1.axis('equal')\n",
    "ax1.set_title(\"Plochy budov v oblasti s teplotou pod 21°C\", fontsize=17, color=colors[1], fontweight=\"bold\", x=0.6)\n",
    "plt.axis(\"equal\")\n",
    "plt.tight_layout()\n",
    "fig1.savefig('09_Pie_below21.svg', transparent=True, bbox_inches='tight')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Area ratio of dominant buildings in the tiles with temperature below 21°C <br>\n",
    "- 46.1% stavba pro rodinnou rekreaci\n",
    "- 23.0% rodinny dum"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(12, 12))\n",
    "ax.set_axis_off()\n",
    "gdf_bh[gdf_bh['zpusob_vyu'] == 'stavba pro rodinnou rekreaci'].centroid.plot(ax=ax, markersize=10, alpha=0.5, color=colors[0])\n",
    "\n",
    "ctx.add_basemap(ax, crs=gdf_bh.crs.to_string(), source=ctx.providers.Stamen.Terrain)\n",
    "plt.tight_layout()\n",
    "fig.savefig('10_Below21_buildings_map.png', bbox_inches='tight', dpi=250)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Location of buildings \"stavba pro rodinnou rekreaci\" is mainly on the outskirts and remote places of Brno."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "buildings_below21_cnt = df_bh_grid[df_bh_grid[\"Temperature\"] < 21 ]['cnt'].sum(axis=0).sort_values(ascending=True)[-18:]\n",
    "\n",
    "fig, ax = plt.subplots(1, 1, figsize=(9, 10))\n",
    "ax.set_xlabel('')\n",
    "color = plt.cm.plasma(np.linspace(0.32,0.10,18))\n",
    "buildings_below21_cnt.plot(kind='barh', ax=ax, zorder=2, color=color)\n",
    "ax.grid(axis='x', linestyle='--', linewidth=\"0.8\", color=\"k\", alpha=0.3, zorder=0)\n",
    "\n",
    "ax.spines[['top','right','bottom', 'left']].set_visible(False)\n",
    "ax.yaxis.set_ticks_position('none')\n",
    "ax.xaxis.set_ticks_position('none')\n",
    "ax.set_ylabel(None)\n",
    "ax.set_yticklabels(ax.get_yticklabels(), fontsize=14)\n",
    "xlabels = ax.get_xticks().astype('i')[1:-1]\n",
    "ax.set_xticks(xlabels)\n",
    "ax.set_xticklabels(xlabels, fontsize=14)\n",
    "ax.tick_params(axis='y', pad=12)\n",
    "\n",
    "for i, v in enumerate(buildings_below21_cnt.values):\n",
    "    plt.text(v+25, i, str(int(v)), color=color[i], va=\"center\", fontsize=14)\n",
    "\n",
    "fig.savefig('11_Barh_below21.svg', bbox_inches='tight')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Absolute count of individual buildings inside the tiles with temperature below 21°C"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Buildings in tiles with temperature between 21°C÷25°C"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_between = pd.Series(df_bh_grid[(df_bh_grid[\"Temperature\"] <= 25) & (df_bh_grid[\"Temperature\"] >= 21)]['b_area_ratio'].sum(axis=0) * 100 / df_bh_grid[(df_bh_grid[\"Temperature\"] <= 25) & (df_bh_grid[\"Temperature\"] >= 21)]['b_area_ratio'].sum(axis=0).sum(axis=0), name=\"Occupy\").to_frame()\n",
    "df_between = df_between[df_between.values>2].sort_values(by='Occupy', ascending=False)\n",
    "df_between"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note: sum of area percentage in \"Occupy\" column is not 100% because buildings with area < 2% have been dropped"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Average height calculated from weighted mean height values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "avg_height_between = gdf_grid_extended.loc[df_bh_grid[(df_bh_grid[\"Temperature\"] <= 25) & (df_bh_grid[\"Temperature\"] >= 21)].index, \"WM_Height\"].sum()/len(df_bh_grid[(df_bh_grid[\"Temperature\"] <= 25) & (df_bh_grid[\"Temperature\"] >= 21)].index)\n",
    "print(\"Average height of tiles with temperature between 21°C÷25°C is %.2f m.\" %avg_height_between)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "labels = df_between.index.values\n",
    "values = df_between.values.flatten()\n",
    "explode = np.size(df_between.values)*(0.01,)\n",
    "explode = list(explode)\n",
    "explode[0] = 0.08\n",
    "explode[1] = 0.08\n",
    "explode[6] = 0.12\n",
    "explode = tuple(explode)\n",
    "\n",
    "colors = plt.cm.plasma(np.linspace(0.5,0.2,np.size(df_between.values)))\n",
    "fig1, ax1 = plt.subplots(figsize=(8,5))\n",
    "patches, labels, pct_texts = ax1.pie(values, labels=labels, labeldistance=1.05, autopct='%1.1f%%', pctdistance=0.8,\n",
    "                                explode=explode, colors=colors, shadow=False, startangle=20, rotatelabels =True, textprops={'fontsize': 12})\n",
    "\n",
    "for label, pct_text in zip(labels, pct_texts):\n",
    "    pct_text.set_rotation(label.get_rotation())\n",
    "\n",
    "for label, pct_text in zip(labels, pct_texts):\n",
    "    label.set_rotation(0)\n",
    "    pct_text.set_color('white')\n",
    "\n",
    "ax1.axis('equal')\n",
    "ax1.set_title(\"Plochy budov v oblasti s teplotou medzi 21°C÷25°C\", fontsize=17, color=colors[1], fontweight=\"bold\", x=0.55)\n",
    "plt.axis(\"equal\")\n",
    "plt.tight_layout()\n",
    "fig1.savefig('12_Pie_between.svg', bbox_inches='tight')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Area ratio of dominant buildings in the tiles with temperature between 21°C÷25°C <br>\n",
    "- 23.6% rodinny dum\n",
    "- 15.7% bytovy dum\n",
    "- 11.5% jina stavba\n",
    "- 10.6% stavba obcanskeho vybaveni"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(12, 12))\n",
    "ax.set_axis_off()\n",
    "gdf_bh[gdf_bh['zpusob_vyu'] == 'rodinný dům'].centroid.plot(ax=ax, markersize=10, alpha=0.4, color=colors[0])\n",
    "gdf_bh[gdf_bh['zpusob_vyu'] == 'bytový dům'].centroid.plot(ax=ax, markersize=10, alpha=0.4, color=colors[1])\n",
    "gdf_bh[gdf_bh['zpusob_vyu'] == 'objekt k bydlení'].centroid.plot(ax=ax, markersize=10, alpha=0.4, color=colors[6])\n",
    "\n",
    "ctx.add_basemap(ax, crs=gdf_bh.crs.to_string(), source=ctx.providers.Stamen.Terrain)\n",
    "plt.tight_layout()\n",
    "fig.savefig('13_Between_buildings_map.png', bbox_inches='tight', dpi=250)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Buildings \"rodinny dum\", \"bytovy dum\" and \"objekt k bydleni\" are located in the center of Brno, as well as on the outskirts."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "buildings_between_cnt = df_bh_grid[(df_bh_grid[\"Temperature\"] <= 25) & (df_bh_grid[\"Temperature\"] >= 21)]['cnt'].sum(axis=0).sort_values(ascending=True)[-17:]\n",
    "\n",
    "fig, ax = plt.subplots(1, 1, figsize=(9, 10))\n",
    "ax.set_xlabel('')\n",
    "color = plt.cm.plasma(np.linspace(0.22,0.5,17))\n",
    "buildings_between_cnt.plot(kind='barh', ax=ax, zorder=2, color=color)\n",
    "ax.grid(axis='x', linestyle='--', linewidth=\"0.8\", color=\"k\", alpha=0.2, zorder=0)\n",
    "\n",
    "ax.spines[['top','right','bottom', 'left']].set_visible(False)\n",
    "ax.yaxis.set_ticks_position('none')\n",
    "ax.xaxis.set_ticks_position('none')\n",
    "ax.set_ylabel(None)\n",
    "ax.set_yticklabels(ax.get_yticklabels(), fontsize=14)\n",
    "xlabels = ax.get_xticks().astype('i')[1:-1]\n",
    "ax.set_xticks(xlabels)\n",
    "ax.set_xticklabels(xlabels, fontsize=14)\n",
    "ax.tick_params(axis='y', pad=12)\n",
    "\n",
    "for i, v in enumerate(buildings_between_cnt.values):\n",
    "    plt.text(v+200, i, str(int(v)), color=color[i], va=\"center\", fontsize=14)\n",
    "\n",
    "fig.savefig('14_Barh_between.svg', bbox_inches='tight')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Absolute count of individual buildings inside the tiles with temperature between 21°C÷25°C"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Buildings in tiles with temperature above 25°C"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_above25 = pd.Series(df_bh_grid[df_bh_grid[\"Temperature\"] > 25 ]['b_area_ratio'].sum(axis=0) * 100 / df_bh_grid[df_bh_grid[\"Temperature\"] > 25 ]['b_area_ratio'].sum(axis=0).sum(axis=0), name=\"Occupy\").to_frame()\n",
    "df_above25 = df_above25[df_above25.values>2].sort_values(by='Occupy', ascending=False)\n",
    "df_above25"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note: sum of area percentage in \"Occupy\" column is not 100% because buildings with area < 2% have been dropped"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Average height calculated from weighted mean height values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "avg_height_above25 = gdf_grid_extended.loc[df_bh_grid[df_bh_grid[\"Temperature\"] > 25].index, \"WM_Height\"].sum() / len(df_bh_grid[df_bh_grid[\"Temperature\"] > 25 ].index)\n",
    "print(\"Average height of tiles with temperature above 25°C is %.2f m.\" %avg_height_above25)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "labels = df_above25.index.values\n",
    "values = df_above25.values.flatten()\n",
    "explode = np.size(df_above25.values)*(0.02,)\n",
    "explode = list(explode)\n",
    "explode[0] = 0.08\n",
    "explode[3] = 0.12\n",
    "explode = tuple(explode)\n",
    "\n",
    "colors = plt.cm.plasma(np.linspace(0.9,0.75,np.size(df_above25.values)))\n",
    "fig1, ax1 = plt.subplots(figsize=(8,5))\n",
    "patches, labels, pct_texts = ax1.pie(values, labels=labels, labeldistance=1.05, autopct='%1.1f%%', pctdistance=0.8,\n",
    "                                explode=explode, colors=colors, shadow=False, startangle=45, rotatelabels =True, textprops={'fontsize': 12})\n",
    "\n",
    "for label, pct_text in zip(labels, pct_texts):\n",
    "    pct_text.set_rotation(label.get_rotation())\n",
    "\n",
    "for label, pct_text in zip(labels, pct_texts):\n",
    "    label.set_rotation(0)\n",
    "\n",
    "ax1.axis('equal')\n",
    "ax1.set_title(\"Plochy budov v oblasti s teplotou nad 25°C\", fontsize=17, color=colors[1], fontweight=\"bold\", x=0.52)\n",
    "plt.axis(\"equal\")\n",
    "plt.tight_layout()\n",
    "fig1.savefig('15_Pie_above25.svg', bbox_inches='tight')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Area ratio of dominant buildings in the tiles with temperature above 25°C <br>\n",
    "- 50.9% stavba pro vyrobu a skladovani\n",
    "- 12.0% stavba obcanskeho vybaveni\n",
    "- 11.5% jina stavba"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(12, 12))\n",
    "ax.set_axis_off()\n",
    "gdf_bh[gdf_bh['zpusob_vyu'] == 'stavba pro výrobu a skladování'].centroid.plot(ax=ax, markersize=15, alpha=0.6, color=colors[0])\n",
    "gdf_bh[gdf_bh['zpusob_vyu'] == 'průmyslový objekt'].centroid.plot(ax=ax, markersize=15, alpha=0.6, color=colors[3])\n",
    "\n",
    "ctx.add_basemap(ax, crs=gdf_bh.crs.to_string(), source=ctx.providers.Stamen.Terrain)\n",
    "plt.tight_layout()\n",
    "fig.savefig('16_Above25_buildings_map.png', bbox_inches='tight', dpi=250)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Industrial buildings \"stavba pro vyrobu a skladovani\" and \"prumyslovy objekt\" are mainly located in the east of the Brno."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "buildings_above25_cnt = df_bh_grid[df_bh_grid[\"Temperature\"] > 25]['cnt'].sum(axis=0).sort_values(ascending=True)[-17:]\n",
    "\n",
    "fig, ax = plt.subplots(1, 1, figsize=(9, 10))\n",
    "ax.set_xlabel('')\n",
    "color = plt.cm.plasma(np.linspace(0.7,0.88,17))\n",
    "buildings_above25_cnt.plot(kind='barh', ax=ax, zorder=2, color=color)\n",
    "ax.grid(axis='x', linestyle='--', linewidth=\"0.8\", color=\"k\", alpha=0.2, zorder=0)\n",
    "\n",
    "ax.spines[['top','right','bottom', 'left']].set_visible(False)\n",
    "ax.yaxis.set_ticks_position('none')\n",
    "ax.xaxis.set_ticks_position('none')\n",
    "ax.set_ylabel(None)\n",
    "ax.set_yticklabels(ax.get_yticklabels(), fontsize=14)\n",
    "xlabels = ax.get_xticks().astype('i')[1:-2]\n",
    "ax.set_xticks(xlabels)\n",
    "ax.set_xticklabels(xlabels, fontsize=14)\n",
    "ax.tick_params(axis='y', pad=12)\n",
    "\n",
    "for i, v in enumerate(buildings_above25_cnt.values):\n",
    "    plt.text(v+10, i, str(int(v)), color=color[i], va=\"center\", fontsize=14)\n",
    "\n",
    "fig.savefig('17_Barh_above25.svg', bbox_inches='tight')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Absolute count of individual buildings inside the tiles with temperature above 25°C"
   ]
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "e404b59586357c814bc0d3940e75d6763c00a48753b225b81f7716971b8e1741"
  },
  "kernelspec": {
   "display_name": "Python 3.8.10 64-bit",
   "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.8.10"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}