############################################################################ # # Name: network_sim.py # Python: 3.3.2 on win32 # Autor: Adrian Haas # ############################################################################ # Simulate a network ############################################################################ # # Simulating a two level star network with backbone nodes. The planning # area is restricted to Zurich. The nodes were genereated randomly. The # linking is done by nearest neighbour search. Design rule: First level # nodes connected to core nodes, second level to first level nodes. The # output is a KML file, which is pushed to map.geo.admin.ch. # ############################################################################ import matplotlib.path as mplPath from scipy import spatial import numpy as np import matplotlib.pyplot as plt import random ############################################################################ #Preparation ############################################################################ #setup node and link tables #definitions: # node[[x_coor, y_corr, type (0 to 2), master, # number of nodes nodes]] # link[[nodea_index,nodeb_index]] node=[] link=[] #define planning area for random node generation #polygons->copy from KML files (map.geo.admin.ch) #inside zurich city poly=np.array([[8.566381817762186,47.34491891234041], [8.58268665804545,47.35231451660678], [8.568695447215644,47.36666476584499], [8.573203153046458,47.379393753843715], [8.552758650817049,47.389667289966916], [8.555690043965356,47.40313164581354], [8.571842856897792,47.402614586627315], [8.595728060078226,47.39212341655013], [8.602051774682427,47.402674228188495], [8.572337449099006,47.41340385653124], [8.563393705792901,47.42914261133091], [8.524278886771073,47.43509334399559], [8.488060185590825,47.42697247082547], [8.505859915300471,47.41565580908032], [8.502500383029895,47.40669179437359], [8.481427804068923,47.41318029307511], [8.46820462330658,47.38721343635927], [8.48950381893336,47.36579334964949], [8.505397611772824,47.352874559506404], [8.5197959976295,47.357957831983704], [8.515654158830875,47.33640811168292], [8.540831616950427,47.338149397154154], [8.533457134444722,47.352791816290406], [8.53555407761833,47.364645548283995], [8.545133089562053,47.36689291420201], [8.552010767233424,47.353873729579696], [8.559368917874803,47.35128387800327], [8.566381817762186,47.34491891234041]]) #outside wood poly_wood=np.array([[8.511355389184171,47.39959967058902], [8.519054211054314,47.40033750583141], [8.529942296782996,47.39830126597944], [8.532350429151934,47.39944785806326], [8.53018806548786,47.40383090186106], [8.522750073810741,47.409432850273305], [8.518283128698785,47.41131863604629], [8.5136888135421,47.41019210594631], [8.512674346149316,47.40916711148981], [8.511858632926515,47.40814026289521], [8.513018273332901,47.40651038850534], [8.508618976624845,47.40187379264983], [8.511355389184171,47.39959967058902]]) #define square around planning area for random coordinate generation a=np.min(poly[:,0]) b=np.max(poly[:,0]) c=np.min(poly[:,1]) d=np.max(poly[:,1]) #polygon objects for inside check bbPath = mplPath.Path(poly) bbPath_wood = mplPath.Path(poly_wood) ############################################################################ #Define core nodes ############################################################################ node.append([8.4885,47.394,0,0,0]) node.append([8.54404,47.41872,0,1,0]) node.append([8.54551,47.37427,0,2,0]) core_node_nr=len(node) core_node=node.copy() ############################################################################ #Add first level nodes to network ############################################################################ i=0 while i<100: x=random.uniform (a,b) y=random.uniform (c,d) #redo random generation until inside planning area while (bbPath.contains_point((x, y))==0 or bbPath_wood.contains_point((x, y))==1): x=random.uniform (a,b) y=random.uniform (c,d) node.append([x,y,1,0,0]) i+=1 #connect to core, use nearest neighbour search corr_table=np.array(core_node)[:,:2] j=0 for row in node: if row[2]==1: pt=[row[0],row[1]] distance,index = spatial.KDTree(corr_table).query(pt) link.append([index,j]) row[3]=index #update master field j+=1 ############################################################################ #Add second level nodes to network ############################################################################ corr_table=np.array(node.copy())[core_node_nr:,:2] i=0 while i<300: x=random.uniform (a,b) y=random.uniform (c,d) #redo random generation until inside planning area while (bbPath.contains_point((x, y))==0 or bbPath_wood.contains_point((x, y))==1): x=random.uniform (a,b) y=random.uniform (c,d) node.append([x,y,2,0,0]) i+=1 #connect to first level nodes, use nearest neighbour search j=0 for row in node: if row[2]==2: pt=[row[0],row[1]] distance,index = spatial.KDTree(corr_table).query(pt) link.append([index+core_node_nr,j]) row[3]=node[index+core_node_nr][3] #define master same as far end j+=1 #network analyis: count connected nodes for each core node stat=np.array(node)[:,3] node[0][4]=len(np.where(stat==0)[0]) node[1][4]=len(np.where(stat==1)[0]) node[2][4]=len(np.where(stat==2)[0]) ############################################################################ #matplotlib plot ############################################################################ #draw city polygon i=0 while i<len(poly)-1: plt.plot([poly[i,0], poly[i+1,0]],[poly[i,1], poly[i+1,1]]) i+=1 plt.plot([poly[i,0], poly[0,0]],[poly[i,1], poly[0,1]]) #draw nodes for row in node: if row[2]==0: plt.plot(row[0],row[1],"ko") elif row[2]==1: plt.plot(row[0],row[1],"r^") else: plt.plot(row[0],row[1],"b^") #draw links for row in link: plt.plot([node[row[0]][0],node[row[1]][0]], [node[row[0]][1],node[row[1]][1]],color='k') plt.show() ############################################################################ #create kml file for view on map.geo.admin.ch ############################################################################ #define kml functions def kopf(): #define styles in kml icon1="http://www.engineering-tools.ch/kml/core_node.png" icon2="http://www.engineering-tools.ch/kml/node1.png" icon3="http://www.engineering-tools.ch/kml/node2.png" s ="""<?xml version='1.0' encoding='UTF-8'?> <kml xmlns='http://www.opengis.net/kml/2.2'> <Document> <Style id="line"> <LineStyle> <color>7f00ffff</color> <width>4</width> </LineStyle> </Style> <Style id="core node"> <LabelStyle> <color>7fff0000</color> <scale>0.8</scale> </LabelStyle> <IconStyle> <Icon> <href>"""+icon1+"""</href></Icon> <scale>1</scale> <heading>0</heading> <hotSpot x='0.4' y='-0.4' xunits='fraction' yunits='fraction'/> </IconStyle> </Style> <Style id="node1"> <LabelStyle> <color>7fff0000</color> <scale>0.8</scale> </LabelStyle> <IconStyle> <Icon> <href>"""+icon2+"""</href> </Icon> <scale>1</scale> <heading>0</heading> <hotSpot x='0.4' y='-0.4' xunits='fraction' yunits='fraction'/> </IconStyle> </Style> <Style id="node2"> <LabelStyle> <color>7fff0000</color> <scale>0.8</scale> </LabelStyle> <IconStyle> <Icon> <href>"""+icon3+"""</href> </Icon> <scale>1</scale> <heading>0</heading> <hotSpot x='0.4' y='-0.4' xunits='fraction' yunits='fraction'/> </IconStyle> </Style> <Placemark> <name>Simuliertes Netzwerk</name> </Placemark> """ return s def kml_node(nr,x,y,typ,core,anz): #add placemark for the node element s=""" <Placemark> <name>"""+str(nr)+"""</name> <description><![CDATA[<h2>Zürich</h2> <br> Typ:"""+str(typ)+"""<br> Core: Node"""+str(core)+"""<br>""" if typ=="core node": s+="Connected Nodes:"""+str(anz)+"""<br>""" s+="""]]> </description> <styleUrl>#"""+str(typ)+"""</styleUrl> <Point> <coordinates>"""+str(x)+""","""+str(y)+""",0</coordinates> </Point> </Placemark> """ return s def kml_link(nr,x1,y1,x2,y2,na,nb): #add placemark for the link element s=""" <Placemark> <name>"""+str(nr)+"""</name> <description><![CDATA[Site A: Node"""+str(na)+"""<br> Site B: Node"""+str(nb)+"""<br>]]> </description> <styleUrl>#line</styleUrl> <LineString> <extrude>1</extrude> <tessellate>1</tessellate> <altitudeMode>absolute</altitudeMode> <coordinates>"""+str(x1)+""", """+str(y1)+""",0 """+str(x2)+""", """+str(y2)+""",0 </coordinates> </LineString> </Placemark> """ return s def ende(): #add the tail s=""" </Document> </kml> """ return s def createfile(s,name): #output outfile=open(name+".kml","w") outfile.write(s) #generate the file translate_type=["core node","node1","node2"] s=kopf() i=0 for row in node: s+=kml_node("Node"+str(i),row[0],row[1],translate_type[row[2]], row[3],row[4]) i+=1 i=0 for row in link: s+=kml_link("Link"+str(i),node[row[0]][0],node[row[0]][1], node[row[1]][0],node[row[1]][1],row[0],row[1]) i+=1 s+=ende() createfile(s,"zh_network")