|
import argparse |
|
import requests |
|
import json |
|
import geopandas as gpd |
|
import pandas as pd |
|
from shapely import affinity |
|
from shapely.geometry import Polygon, Point, MultiPolygon |
|
from matplotlib import pyplot as plt |
|
import numpy as np |
|
from osmtogeojson import osmtogeojson |
|
|
|
|
|
#0) Get user input |
|
cities = ["Zurich", "New_York_City","New_Delhi", "Berlin", "Barcelona", "Tokyo"] |
|
|
|
citiesBbox_wgs = { |
|
"Zurich": "(47.362414,8.513328,47.385392,8.559779)", #zurich |
|
"New_York_City": "(40.713196, -74.017528, 40.739165, -73.967435)",#nyc |
|
"New_Delhi": "(28.599522,77.194575,28.62432, 77.240522)", #new delhi |
|
"Berlin":"(52.503223,13.350635,52.547340,13.447870)", #berlin |
|
"Barcelona": "(41.375642,2.129504,41.401074,2.180285)", #barcelona |
|
"Tokyo": "(35.640169,139.712369,35.704226,139.801369)" #tokyo |
|
} |
|
|
|
parser = argparse.ArgumentParser() |
|
parser.add_argument("-c","--city", required=True, help="Select a city to organize", choices=cities) |
|
city = parser.parse_args() |
|
userinput = vars(parser.parse_args())["city"] |
|
|
|
#1) Get data from OSM |
|
def getData(citycode): |
|
overpass_query = """ |
|
[out:json][timeout:25]; |
|
(way["building"]"""+ citycode +"""; |
|
relation["building"]"""+ citycode +"""; |
|
); |
|
out body; |
|
>; |
|
out skel qt; |
|
""" |
|
overpass_url = "http://overpass-api.de/api/interpreter?" |
|
response = requests.get(overpass_url, params={'data': overpass_query}) |
|
print(response) |
|
response.raise_for_status() |
|
return response.json() |
|
|
|
|
|
def getProcessedData(citycode): |
|
try: |
|
r = getData(citycode) |
|
return r |
|
except requests.exceptions.ConnectionError: |
|
print("Check your internet connection") |
|
except requests.exceptions.HTTPError: |
|
print("Sorry, there was a HTTPError in the HTTP request returned") |
|
|
|
|
|
#2) Sort and prepare geodataframe by adding necessary values for the geometric translation |
|
def prepareDf(data, sortFeature): |
|
city = [] |
|
for i in data["features"]: |
|
if "coordinates" in i["geometry"] and i["geometry"]["type"] == "Polygon": |
|
city.append(MultiPolygon([Polygon(poly) for poly in i["geometry"]["coordinates"]])) |
|
|
|
city = gpd.GeoSeries(city) |
|
city.crs = {'init': 'epsg:4326', 'no_defs': True} |
|
city_proj = city.to_crs({'init': 'epsg:3857'}) |
|
|
|
df = gpd.GeoDataFrame(pd.concat([city_proj.area, city_proj.length], axis=1), geometry = city_proj) |
|
df.rename(columns = {0: "area", 1: "length"}, inplace=True) |
|
df = df[df[sortFeature] != 0.0].sort_values(sortFeature).reset_index(drop=True) |
|
df["lengthX"] = (df.bounds["maxx"] - df.bounds["minx"]).abs() |
|
df["lengthY"] = (df.bounds["maxy"] - df.bounds["miny"]).abs() |
|
df["cumulative_X"] = df["lengthX"].cumsum() #cumulative length |
|
l_guess = max(df["cumulative_X"])*0.01 |
|
return df, l_guess |
|
|
|
#3) find optimized total width and height (ideal width / height ratio around 0.7) |
|
def findWidth(df, l_guess, counter, max_iter): |
|
tolerance = 0.01 #tolerated offset to expected ratio |
|
r_expected = 0.7 #expected ratio |
|
h = 0 |
|
counter += 1 |
|
df = df.copy() |
|
df["row"] = [np.floor(i/l_guess) for i in df["cumulative_X"]] |
|
|
|
for k in range(int(max(df["row"])+1)): |
|
vals = [j for j,i in enumerate(df["row"]) if i == k] |
|
h += max(df['lengthY'][vals])*1.1 |
|
df.loc[vals, "new_X"] = df['lengthX'][vals].cumsum().shift(1, fill_value=0) |
|
df.loc[vals, "new_Y"] = h |
|
|
|
r = l_guess/h #calculated ratio |
|
|
|
if abs(r_expected - r) <= tolerance or counter == max_iter: |
|
print("optimal length is:", l_guess, "height is:", h, "ratio:", r) |
|
return df |
|
elif r < r_expected: #l_guess too small |
|
l_guess += l_guess * 0.15 |
|
return findWidth(df,l_guess,counter,max_iter) |
|
else: #l_guess too big |
|
l_guess -= l_guess * 0.15 |
|
return findWidth(df,l_guess,counter, max_iter) |
|
|
|
#4) Create translation vector and shift geometries |
|
def translate(df): |
|
df["translateY"] = -(df.bounds["miny"]) + df["new_Y"] |
|
df["translateX"] = -df.bounds["minx"] + df["new_X"] |
|
df["geometry"] = df.apply(lambda x: affinity.translate(x.geometry, xoff = x.translateX, yoff=x.translateY), axis=1) |
|
return df |
|
|
|
#4) Plot reorganized map and save as image |
|
def draw_city(gdfFinal): |
|
f, ax = plt.subplots(1, figsize = (10,15)) |
|
ax = gdfFinal.plot(ax = ax, color='gray', linewidth =0.2) |
|
lims = plt.axis("equal") |
|
plt.axis('off') |
|
plt.savefig(userinput + ".png") |
|
plt.show() |
|
|
|
def main(): |
|
data = osmtogeojson.process_osm_json(getProcessedData(citiesBbox_wgs[userinput])) |
|
data_prep = prepareDf(data, "length") |
|
data_trans = findWidth(data_prep[0], data_prep[1], 0, 10) |
|
draw_city(translate(data_trans)) |
|
|
|
main() |
