1+ import matplotlib .pyplot as plt
2+ import networkx as nx
3+ import xml .etree .ElementTree as ET
4+ import geopandas as gpd
5+ import numpy as np
6+ from shapely .ops import unary_union
7+ from shapely .geometry import Polygon
8+ import math
9+
10+ # Load data from XML
11+ tree = ET .parse ('Tools\gis2graph\graph_files\King_county_NG911.graphml' )
12+ root = tree .getroot ()
13+
14+ # Initialize the network graph
15+ G = nx .DiGraph ()
16+
17+ # Type attribute key
18+ type_attribute_key = 'type'
19+ # Namespace from root
20+ nsmap = {'xmlns' : 'http://graphml.graphdrawing.org/xmlns' }
21+
22+ node_positions = {}
23+ psap_ids = []
24+ psap_coordinates = []
25+ calr_nodes = []
26+ caller_region_patches = []
27+
28+ region_string = [] # stores the strings of squares that will be set as caller region attributes
29+
30+ out_file_name = "King_county_NG911"
31+
32+ square_counter = 0
33+
34+ psap_layer = gpd .read_file ("Tools\gis2graph\GIS_data\Layers\PSAP_layer.gpkg" )
35+ provisioning_layer = gpd .read_file ("Tools\gis2graph\GIS_data\Layers\Provisioning_layer.gpkg" )
36+
37+ # create series of boolean values denoting whether geometry is within King County
38+ psap_within_kc = psap_layer .within (provisioning_layer .iloc [21 ].geometry )
39+
40+ # create new GeoDataFrames of just items located within King County using series above
41+ kc_psap = psap_layer .drop (np .where (psap_within_kc == False )[0 ])
42+
43+
44+ # area_multiplier gets multiplied by the smallest psap area to determine size of the squares in the grid
45+ area_multiplier = 10
46+
47+ # Some of the data from the state had messed up psap names. This line fixes them so they can be consolidated
48+ kc_psap .loc [[265 , 262 , 185 ], 'DsplayName' ] = "King County Sheriff's Office - Marine Patrol"
49+
50+ names = [] # 2 empty lists for storing names and polygons for merging psaps together
51+ polys = []
52+ es_nguid = [] # list for storing the nguids
53+
54+ # Loops through and finds all unique names, as well as sorting all the polygons that make up those regions
55+ for n in range (kc_psap .shape [0 ]):
56+ if (kc_psap .iloc [n ].DsplayName ) in names :
57+ polys [names .index (kc_psap .iloc [n ].DsplayName )].append (kc_psap .iloc [n ].geometry )
58+ else :
59+ names .append (kc_psap .iloc [n ].DsplayName )
60+ polys .append ([kc_psap .iloc [n ].geometry ])
61+ es_nguid .append (kc_psap .iloc [n ].ES_NGUID )
62+
63+ # Takes the lists of polygons, and merges them into new polygons for the creation of the merged_kc_psap GeoDataFrame
64+ merged_polys = []
65+ for m in range (len (polys )):
66+ merged_polys .append (unary_union (polys [m ]))
67+
68+ # Create a new GeoDataFrame with the unique names and merged geometries
69+ merged_kc_psap = gpd .GeoDataFrame ({'DisplayName' : names , 'geometry' : merged_polys , 'ES_NGUID' : es_nguid }, crs = kc_psap .crs )
70+
71+ # Find the area of the smallest merged psap, use that to determine the square size
72+ areas = merged_kc_psap .area
73+ side_length = math .sqrt (areas .min () * area_multiplier )
74+
75+ # Creates a grid of squares based on the bounds of merged_kc_psaps, and the side_length
76+ xmin , ymin , xmax , ymax = merged_kc_psap .total_bounds
77+ cols = list (np .arange (xmin , xmax + side_length , side_length ))
78+ rows = list (np .arange (ymin , ymax + side_length , side_length ))
79+ squares = []
80+ for x in cols [:- 1 ]:
81+ for y in rows [:- 1 ]:
82+ squares .append (
83+ Polygon ([(x , y ), (x + side_length , y ), (x + side_length , y + side_length ), (x , y + side_length )]))
84+ grid = gpd .GeoDataFrame ({'geometry' : squares }, crs = kc_psap .crs )
85+ kc_psap .plot ()
86+
87+ #this code iterates through graph ML and finds nodes + adds thenm to graph
88+ for node in root .findall ('.//{http://graphml.graphdrawing.org/xmlns}node' ):
89+ type_element = node .find (f'.//{{{ nsmap ["xmlns" ]} { type_attribute_key } )
90+ if type_element is not None and type_element .text == 'PSAP' :
91+ node_id = node .get ('id' )
92+ psap_ids .append (node_id )
93+
94+ for data in node .findall (f'.//{{{ nsmap ["xmlns" ]} , namespaces = nsmap ):
95+ if data .get ('key' ) == 'x' :
96+ node_x = float (data .text )
97+ psap_coordinates .append (node_x )
98+ elif data .get ('key' ) == 'y' :
99+ node_y = float (data .text )
100+ psap_coordinates .append (node_y )
101+
102+ if 'node_x' in locals () and 'node_y' in locals ():
103+ G .add_node (node_id )
104+ node_positions [node_id ] = (node_x , node_y )
105+ G .nodes [node_id ]['pos' ] = (node_x , node_y )
106+ G .nodes [node_id ]['color' ] = 'cyan'
107+ #change EMS with FIRE or LAW to see fire/police nodes
108+ elif type_element is not None and type_element .text == 'EMS' :
109+ node_id = node .get ('id' )
110+
111+ for data in node .findall (f'.//{{{ nsmap ["xmlns" ]} , namespaces = nsmap ):
112+ if data .get ('key' ) == 'x' :
113+ node_x = float (data .text )
114+ if data .get ('key' ) == 'y' :
115+ node_y = float (data .text )
116+
117+ if 'node_x' in locals () and 'node_y' in locals ():
118+ G .add_node (node_id )
119+ node_positions [node_id ] = (node_x , node_y )
120+ G .nodes [node_id ]['pos' ] = (node_x , node_y )
121+ G .nodes [node_id ]['color' ] = 'blue'
122+
123+
124+
125+ #find the segments of caller regions and find average center to mark as point
126+ for node in root .findall ('.//{http://graphml.graphdrawing.org/xmlns}node' ):
127+ type_element = node .find (f'.//{{{ nsmap ["xmlns" ]} { type_attribute_key } )
128+ if type_element is not None and type_element .text == 'CALR' :
129+ node_id = node .get ('id' )
130+ calr_nodes .append (node_id )
131+
132+ segments_data = node .find (f'.//{{{ nsmap ["xmlns" ]} )
133+ if segments_data is not None :
134+ segments_str = segments_data .text
135+ segments_str = segments_str .replace ("[(" , "" ).replace (")]" , "" )
136+ segments_list = segments_str .split ("), (" )
137+
138+ # Convert segments into a list of coordinate tuples
139+ coordinates = [tuple (map (float , segment .split (", " ))) for segment in segments_list ]
140+
141+ # Calculate the average position as the node position- this should be the center of the region- can map
142+ # onto the graph
143+ avg_x = sum (coord [0 ] for coord in coordinates ) / len (coordinates )
144+ avg_y = sum (coord [1 ] for coord in coordinates ) / len (coordinates )
145+
146+ G .add_node (node_id )
147+ node_positions [node_id ] = (avg_x , avg_y )
148+ G .nodes [node_id ]['pos' ] = (avg_x , avg_y )
149+ G .nodes [node_id ]['color' ] = 'green'
150+
151+
152+
153+ for edge in root .findall ('.//{http://graphml.graphdrawing.org/xmlns}edge' ):
154+ edge_source = edge .get ('source' )
155+ edge_target = edge .get ('target' )
156+
157+ # Check if both source and target nodes are in G
158+ if edge_source in G .nodes () and edge_target in G .nodes ():
159+ G .add_edge (edge_source , edge_target )
160+ # Create the combined plot
161+ fig , ax = plt .subplots (figsize = (10 , 10 ))
162+
163+ # Plot the network graph
164+ edge_alpha = 0.4
165+ default_node_color = 'gray'
166+ node_colors = [G .nodes [node_id ].get ('color' , default_node_color ) for node_id in G .nodes ()]
167+ nx .draw (G , node_size = 30 , pos = node_positions , node_color = node_colors , with_labels = False ,
168+ alpha = edge_alpha , width = 0.3 , ax = ax )
169+
170+ label_pos = {k : (v [0 ], v [1 ] + 0.01 ) for k , v in node_positions .items ()}
171+ labels = {node : node for node in G .nodes ()}
172+ nx .draw_networkx_labels (G , label_pos , labels = labels , font_size = 8 , ax = ax )
173+
174+ kc_psap .plot (ax = ax , color = 'none' , edgecolor = 'black' )
175+ plt .show ()
0 commit comments