#!/usr/bin/env python # coding: utf-8 # In[48]: import geopandas as gpd from pathlib import Path import pandas as pd import matplotlib.pyplot as plt import folium import seaborn as sns import datetime as dt import numpy as np from branca.colormap import linear # In[2]: city = 'Luanda' # In[10]: # read catchment AOI cat = gpd.read_file('AOI/luanda_catchment_level4.shp').to_crs(epsg = 4326) aoi = gpd.read_file('AOI/luanda.shp').to_crs(epsg = 4326) centroid = cat.centroid # In[4]: output_folder = Path('output') # ### Dams and reservoirs # In[5]: dams = gpd.read_file(output_folder / f'{city.lower()}_dams.shp') reservoirs = gpd.read_file(output_folder / f'{city.lower()}_reservoirs.shp') # In[6]: dams['DAM_ID'] = dams['DAM_ID'].astype('int').astype('string') # In[9]: # Create a folium map centered on the AOI mymap = folium.Map(tiles="Cartodb Positron", control_scale = True) bounds = cat.total_bounds.tolist() mymap.fit_bounds([bounds[:2][::-1], bounds[2:][::-1]]) folium.GeoJson( cat, style_function=lambda feature: { 'fillColor': 'transparent', 'color': 'gray', 'weight': 3 },).add_to(mymap) res_tooltip = folium.GeoJsonTooltip(fields=['RES_NAME'], aliases=['Reservoir Name'], localize=True, style=('background-color: white; color: #333333;')) dam_tooltip = folium.GeoJsonTooltip(fields=['DAM_ID'], aliases=['Dam ID'], localize=True, style=('background-color: white; color: #333333;')) # Add blue polygons for reservoirs folium.GeoJson( reservoirs, style_function=lambda feature: { 'fillColor': 'blue', 'color': 'blue', 'weight': 1, }, highlight_function=lambda x: {"fillOpacity": 0.7}, tooltip=res_tooltip, zoom_on_click=True ).add_to(mymap) # Add CircleMarkers for dams folium.GeoJson( dams, name="Dams", marker=folium.CircleMarker(radius=8, fill_color="brown", fill_opacity=0.4, color="brown", weight=1), tooltip = dam_tooltip, highlight_function=lambda x: {"fillOpacity": 0.7}, zoom_on_click=True, ).add_to(mymap) # Add legend manually legend_html = '''
  Legend
  Reservoirs  
  Dams  
''' mymap.get_root().html.add_child(folium.Element(legend_html)) mymap # ### Precipitation # In[20]: df = pd.read_csv(output_folder / 'precip.csv') df['time'] = pd.to_datetime(df['time']) # In[45]: # plot # Create the line chart using Seaborn with connected lines plt.figure(figsize=(12, 6)) sns.lineplot(data=df, x='time', y='precip', linewidth=1, color = '#2b8cbe') # sns.lineplot(data=df, x='time', y=df.precip.mean(), linewidth=1, color = 'grey') # Add labels and title plt.xlabel('') plt.ylabel('') plt.title('Monthly precipitation 1993-2022 (mm)') # plt.xticks(rotation=45) plt.grid(False) plt.show() # ### Runoff # In[46]: var_list = {'ncrun': 'Net cell runoff', 'qs': 'Fast surface and fast subsurface runoff'} stat_list = {'mean': 'mean', 'std': 'standard deviation'} year_list = range(1987, 2017) # In[66]: # Create a folium map centered on the AOI mymap = folium.Map(tiles="Cartodb Positron", control_scale = True) aoi_bounds = cat.total_bounds.tolist() mymap.fit_bounds([aoi_bounds[:2][::-1], aoi_bounds[2:][::-1]]) folium.GeoJson( cat, style_function=lambda feature: { 'fillColor': 'transparent', 'color': 'gray', 'weight': 3 }, control = False).add_to(mymap) # In[67]: for var in var_list: # set colormap if var == 'ncrun': colormap = linear.GnBu_09.scale(json['raster_val'].min(), json['raster_val'].max()) colormap.caption = f"{var_list[var]} (mean / standard deviation)" colormap.add_to(mymap) elif var == 'qs': colormap = linear.YlGnBu_09.scale(json['raster_val'].min(), json['raster_val'].max()) colormap.caption = f"{var_list[var]} (mean / standard deviation)" colormap.add_to(mymap) for stat in stat_list: # read json json = gpd.read_file(output_folder / f'{city.lower()}_{var}_{stat}.geojson').to_crs('epsg:4326') # determine whether to show layer by default show = False if var == 'ncrun' and stat == 'mean': show = True # Add GeoJson layer with choropleth style folium.GeoJson( json, name = f'{var_list[var]} ({stat_list[stat]})', # overlay = False, style_function=lambda feature: { 'fillColor': 'transparent' if pd.isna(feature['properties']['raster_val']) else colormap(feature['properties']['raster_val']), # 'fillColor': colormap(feature['properties']['raster_val']), 'fillOpacity': 0.7, 'weight': 0, }, highlight_function=lambda x: {'fillOpacity': 0.9}, tooltip=folium.GeoJsonTooltip(fields=['raster_val'], aliases = [f'{var_list[var]} ({stat_list[stat]})'], labels=True, sticky=False), show = show ).add_to(mymap) # Add Layer Control to toggle the choropleth layer folium.LayerControl().add_to(mymap) mymap # ### Potential evapotranspiration # In[68]: df = pd.read_csv(output_folder / 'potential_evaporation.csv') df['time'] = pd.to_datetime(df['time']) # In[71]: # plot # Create the line chart using Seaborn with connected lines plt.figure(figsize=(12, 6)) sns.lineplot(data=df, x='time', y='Ep', linewidth=1, color = '#c51b8a') # Add labels and title plt.xlabel('') plt.ylabel('') plt.title('Monthly potential evaporation 1993-2022 (mm/day)') # plt.xticks(rotation=45) plt.grid(False) plt.show() # ### Water pollution # In[72]: ec = pd.read_csv(output_folder / f'{city.lower()}_electrical_conductivity.csv') noxn = pd.read_csv(output_folder / f'{city.lower()}_nitrate_nitrite.csv') bod = pd.read_csv(output_folder / f'{city.lower()}_biological_oxygen_demand.csv') # In[73]: ec['day'] = 1 ec['time'] = pd.to_datetime(ec[['year', 'month', 'day']]) noxn['day'] = 1 noxn['time'] = pd.to_datetime(noxn[['year', 'month', 'day']]) bod['day'] = 1 bod['time'] = pd.to_datetime(bod[['year', 'month', 'day']]) # In[75]: # plot # Create the line chart using Seaborn with connected lines plt.figure(figsize=(12, 6)) sns.lineplot(data=ec, x='time', y='preds', linewidth=1, color = '#e34a33') # Add labels and title plt.xlabel('') plt.ylabel('') plt.title('Monthly electrical conductivity in surface water 1993-2022 (uS/cm)') # plt.xticks(rotation=45) plt.grid(False) plt.show() # In[77]: # plot # Create the line chart using Seaborn with connected lines plt.figure(figsize=(12, 6)) sns.lineplot(data=noxn, x='time', y='preds', linewidth=1, color = '#016c59') # Add labels and title plt.xlabel('') plt.ylabel('') plt.title('Monthly nitrate-nitrite in surface water 1993-2022 (mg/l)') # plt.xticks(rotation=45) plt.grid(False) plt.show() # In[78]: # plot # Create the line chart using Seaborn with connected lines plt.figure(figsize=(12, 6)) sns.lineplot(data=bod, x='time', y='preds', linewidth=1, color = '#1c9099') # Add labels and title plt.xlabel('') plt.ylabel('') plt.title('Monthly biological oxygen demand in surface water 1993-2022 (mg/l)') # plt.xticks(rotation=45) plt.grid(False) plt.show() # ### Residential overlay with livestock # In[79]: glw = gpd.read_file(output_folder / f'{city.lower()}_livestock.geojson').to_crs('epsg:4326') # In[84]: # Create a folium map centered on the AOI mymap = folium.Map(tiles="Cartodb Positron", control_scale = True) aoi_bounds = aoi.total_bounds.tolist() mymap.fit_bounds([aoi_bounds[:2][::-1], aoi_bounds[2:][::-1]]) folium.GeoJson( aoi, style_function=lambda feature: { 'fillColor': 'transparent', 'color': 'gray', 'weight': 3 }, control = False).add_to(mymap) # **Could you please add the built-up area underlay to this map??** # In[85]: # Create a colormap based on the 'YlOrRd' scale and the 'raster_val' column colormap = linear.YlOrRd_09.scale(glw['raster_val'].min(), glw['raster_val'].max()) colormap.caption = "Number of livestocks" colormap.add_to(mymap) # Add GeoJson layer with choropleth style folium.GeoJson( glw, name='Number of livestocks', style_function=lambda feature: { 'fillColor': colormap(feature['properties']['raster_val']), 'fillOpacity': 0.7, 'weight': 0, }, highlight_function=lambda x: {'fillOpacity': 0.9}, tooltip=folium.GeoJsonTooltip(fields=['raster_val'], aliases = ['Number of livestocks'], labels=True, sticky=False), ).add_to(mymap) # Add Layer Control to toggle the choropleth layer folium.LayerControl().add_to(mymap) mymap