<jupyter_start><jupyter_text>For our Sierra Leone example, we have flooding data provided by the GFDRR team (Fathom maps). Here, we match this data on to our OSM road network, and adjust travel speeds accordingly.<jupyter_code>import os, sys
import pandas as pd
import geopandas as gpd
import networkx as nx
from shapely.geometry import Point, MultiPoint
from shapely.wkt import loads
from scipy import spatial
from functools import partial
import pyproj
from shapely.ops import transform
sys.path.append(r'C:\Users\charl\Documents\GitHub\GOST_PublicGoods\GOSTNets\GOSTNets')
import GOSTnet as gn
import rasterio
import numpy as np<jupyter_output>networkx version: 2.3 
osmnx version: 0.9<jupyter_text>Here, we choose the flood return period, as well as set the file locations to the relevant rasters<jupyter_code>flood_return_period = 20

pluv_loc = r'C:\Users\charl\Documents\GOST\SierraLeone\Flooding\SL_pluvial_undefended'

pluvial_raster = os.path.join(pluv_loc, r'SL-PU-{}-1.tif'.format(flood_return_period))

fluv_loc = r'C:\Users\charl\Documents\GOST\SierraLeone\Flooding\SL_fluvial_undefended'

fluvial_raster = os.path.join(fluv_loc, r'SL-FU-{}-1.tif'.format(flood_return_period))<jupyter_output><empty_output><jupyter_text>We read in the raster to memory using rasterio<jupyter_code>fluv = rasterio.open(fluvial_raster)

pluv = rasterio.open(pluvial_raster)<jupyter_output><empty_output><jupyter_text>We create two objects of the raw flood depth arrays, and then find the pair-wise maximum of the two arrays to generate a new array, 'max_depth', which describes the maximum expected depth of flood water, from either pluvial or fluvial floods, that you might expect to occur in DRC over a twenty year time span (where twenty years is the 'return period' of the flooding event)<jupyter_code>fluv_array = fluv.read(1)
pluv_array = pluv.read(1)

max_depth = np.maximum(fluv_array, pluv_array)<jupyter_output><empty_output><jupyter_text>We write this combined raster out to a new raster file, to keep for use later<jupyter_code>out_path = r'C:\Users\charl\Documents\GOST\SierraLeone\Flooding\combined'

final_flood_path = os.path.join(out_path,'combo_%s.tif' % flood_return_period)

meta = pluv.meta.copy()

with rasterio.open(final_flood_path, 'w', **meta) as out:

    out.write_band(1, max_depth)<jupyter_output><empty_output><jupyter_text>Now, we intersect this flood raster with our road network. First we load in the road network:<jupyter_code>net_pth = r'C:\Users\charl\Documents\GOST\SierraLeone\RoadNet'

G = nx.read_gpickle(os.path.join(net_pth, 'final_G.pickle'))<jupyter_output><empty_output><jupyter_text>Then, we intersect the graph's nodes with the combined flood layer to match on the values:<jupyter_code>G_flood = gn.sample_raster(G, final_flood_path, property_name = 'flood_depth')<jupyter_output><empty_output><jupyter_text>We want to go through and make sure no 'bad values' were matched on. These are values above 10m or below 0m. We set these back to 0m (no disruption). These values occur if nodes lie either in permanent water bodies or out of the range of the raster in the no data areas.<jupyter_code>for u, data in G_flood.nodes(data = True):
    fl_val = data['flood_depth']
    if 0 < fl_val < 10:
        pass
    else:
        data['flood_depth'] = 0<jupyter_output><empty_output><jupyter_text>Although a function exists in GOSTnets for breaking a graph for a given node property at a given threshold, in this exercise, we want to prevent travel after depths of 0.5m, but linearly increase travel time for flood depths below this. As such, we will copy and paste the 'disrupt network' function, but make a few local changes for our use case:<jupyter_code>def disrupt_network(G, attr, thresh, fail_value = 99999999):
    
    #### Function for disrupting a graph ####
    # REQUIRED: G - a graph containing one or more nodes and one or more edges
    #           property - the element in the data dictionary for the edges to test
    #           thresh - values of data[property] above this value are disrupted
    #           fail_value - the data['time'] property is set to this value to simulate the removal of the edge
    # RETURNS:  a modified graph with the edited 'time' attribute
    # -------------------------------------------------------------------------#

    G_copy = G.copy()     # copy the input graph

    broken_nodes, disrupted_nodes = [], []    # generate some empty list objects
    
    disruption_dict = {}   # generate an empty dictionary

    for u, data in G_copy.nodes(data = True):   # Iterate through all nodes, one at a time

        if data[attr] > thresh:    # if larger than the threshold...

            broken_nodes.append(u)    # record this as a broken node
            
            disruption_dict[u] = data[attr]  # add node ID and corresponding flood depth to a dictionary
            
        elif 0 < data[attr] < thresh:    # if flood depth is bigger than 0 - but smaller than thresh:
            
            disrupted_nodes.append(u)   # record this as a disrupted node
            
            disruption_dict[u] = data[attr]   # add node ID and corresponding flood depth to a dictionary
        
        else: # If this node is untouched by the flood....
            
            disruption_dict[u] = 0    # ....set its dictionary value to 0
            
    i = 0
    
    for u, v, data in G_copy.edges(data = True):     # starting a second iteration, this time of edges:

        if u in broken_nodes or v in broken_nodes:    # if start or end node is a broken node

            data['time_adj'] = fail_value    # prevent travel along this edge - set travel time to high value
            i+=1
        
        elif u in disrupted_nodes or v in disrupted_nodes: # if start or end node is a disrupted node
            
            depth = max(disruption_dict[u], disruption_dict[v]) # flood depth is max of depth of flood at u and v
            
            orig_time = data['time_adj'] 
            
            data['time_adj'] = orig_time * (1/ ((thresh - depth) / thresh))  # linearly increase travel times 

    print('edges disrupted: %s' % i)
    
    return G_copy<jupyter_output><empty_output><jupyter_text>With our custom edge disruption function written, we can now apply this to the flood graph:<jupyter_code>G_flood = disrupt_network(G_flood, 'flood_depth', 0.5)<jupyter_output>edges disrupted: 7208<jupyter_text>At this point, we re-save our graph, knowing that travel will be disrupted in some cases in this model. In order to compare the disrupted to the non-disrupted network, we will have to generate to OD-matrices.<jupyter_code>gn.save(G_flood, 'final_G_flood', net_pth, , nodes = False,  edges = False)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.08912946759602551, 'fraction_numerical': 0.006228681595728904, 'mean_word_length': 3.4870259481037924, 'pattern_counts': {'":': 0, '<': 47, '<?xml version=': 0, '>': 44, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 1, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '11', 'dir': 'jupyter-structured-clean-dedup', 'id': '332063', 'license': 'MIT', 'n_tokens_mistral': 2172, 'n_tokens_neox': 2062, 'n_words': 731, 'path': 'Implementations/FY20/ACC_SierraLeone - GOSTnets/Step 4.a - Generate Flooding Version of Graph.ipynb', 'repo_name': 'worldbank/GOST_PublicGoods', 'size': '6743'}