tools/simplify_country.py

   1 import argparse
   2 from collections import OrderedDict
   3 import fiona
   4 import shapely.geometry
   5 from fiona.crs import CRS
   6 from fiona.transform import transform_geom
   7 from shapely.wkt import dumps
   8
   9 # data from https://www.naturalearthdata.com/downloads/50m-cultural-vectors/
  10 ZIP_FILE = '/home/philipp/daten/GeoData/naturalearth/v5.1/ne_50m_admin_0_countries.zip'
  11 LAYER_NAME = 'ne_50m_admin_0_countries'
  12
  13
  14 def load_country(country_name: str):
  15     with fiona.open(f'zip://{ZIP_FILE}', layer=LAYER_NAME) as src:
  16         for feature in src:
  17             properties = feature['properties']
  18             if properties['NAME'] == country_name:
  19                 return feature['geometry'], src.crs
  20
  21
  22 def simplify_inner(country, crs):
  23     metric_crs = CRS.from_epsg(32632)  # Western Austria (UTM32)
  24     country = transform_geom(crs, metric_crs, country)
  25     country = shapely.geometry.shape(country)
  26     country = country.buffer(-5000)
  27     country = country.simplify(tolerance=5000)
  28     country = shapely.geometry.mapping(country)
  29     country = transform_geom(metric_crs, crs, country)
  30     return country
  31
  32
  33 def save_country(country, crs, name, shapefile):
  34     schema = {
  35         'geometry': 'Polygon',
  36         'properties': OrderedDict([('name', 'str')])
  37     }
  38     with fiona.open(shapefile, 'w', crs=crs, schema=schema, driver='Shapefile') as c:
  39         feature = {
  40             'geometry': country,
  41             'properties': OrderedDict([('name', name)]),
  42         }
  43         c.write(feature)
  44
  45
  46 def print_wkt(country):
  47     country = shapely.geometry.shape(country)
  48     print(dumps(country, rounding_precision=3))
  49
  50
  51 def main(country_name: str, shapefile: str | None):
  52     country, crs = load_country(country_name)
  53     simplified_country = simplify_inner(country, crs)
  54     if shapefile:
  55         save_country(simplified_country, crs, country_name, shapefile)
  56     print_wkt(simplified_country)
  57
  58
  59 if __name__ == '__main__':
  60     parser = argparse.ArgumentParser(description='Creates a simplified version of a country')
  61     parser.add_argument('--country', default='Austria', help='Country to simplify (e.g. Switzerland)')
  62     parser.add_argument('--shapefile', help='Write result to shape file with this name')
  63     args = parser.parse_args()
  64     main(args.country, args.shapefile)