anisotropi4 · April 20, 2021 21:14
diff --git a/.blocks b/.blocks
 Released license:MIT
 height:780
 border:no

 Towns: Location of urbans areas > 10k population in Great Britain

 Digital boundaries and reference map data:
 Source: Office for National Statistics licensed under the Open Government Licence v.3.0. 
 Source: National Records Scotland data © Crown copyright and database right 2020.
 Contains OS data © Crown copyright and database right 2020.
diff --git a/README.md b/README.md
diff --git a/index.html b/index.html
 <!DOCTYPE html>
 <html>
    <head>
        <title>Towns: Location in GB</title>
        <meta charset="utf-8"/>
        <meta name="viewport" content="width=device-width, initial-scale=1.0">
        <meta http-equiv="Cache-control" content="No-Cache">
        <link rel="shortcut icon" type="image/x-icon" href="docs/images/favicon.ico"/>
        <link rel="stylesheet" href="https://unpkg.com/[email protected]/dist/leaflet.css"
              integrity="sha512-xwE/Az9zrjBIphAcBb3F6JVqxf46+CDLwfLMHloNu6KEQCAWi6HcDUbeOfBIptF7tcCzusKFjFw2yuvEpDL9wQ==" crossorigin=""/>
        <script src="https://unpkg.com/[email protected]/dist/leaflet.js"
                integrity="sha512-gZwIG9x3wUXg2hdXF6+rVkLF/0Vi9U8D2Ntg4Ga5I5BZpVkVxlJWbSQtXPSiUTtC0TjtGOmxa1AJPuV0CPthew=="
                crossorigin=""></script>
        <script src="http://d3js.org/d3.v5.min.js"></script>
        <div id="map" style="width: 1340px; height: 780px"></div>
    </head>
    <body>
        <div id='map'></div>
        <script type="text/javascript">
         var radius = 3;
         var weight = 1;
         var linewidth = 2;
         var log2 = Math.log(2.0);
         var minZoom = 3;
         var maxZoom = 18;
         var map = L.map('map').setView([54.533, -2.53], 6);

         L.tileLayer("http://{s}.tile.openstreetmap.org/{z}/{x}/{y}.png", {
             maxZoom: maxZoom,
             minZoom: minZoom,
             attribution: 'Map data &copy; <a href="https://www.openstreetmap.org/">OpenStreetMap</a> contributors, ' + '<a href="https://creativecommons.org/licenses/by-sa/2.0/">CC-BY-SA</a>, ' + 'Imagery © <a href="https://www.mapbox.com/">Mapbox</a>',
             id: 'osm.standard'
         }).addTo(map);

         d3.json("town_boundaries.geojson").then(function(d) {
             function onEachFeature(feature, layer) {
                 var this_feature = feature.properties;
                 var popupContent;
                 if (this_feature.Town)
                     popupContent = this_feature.Town;
                 else
                     popupContent = this_feature.railway;
                 var lookup = {
                     "bua_code": "BUA code",
                     "bua_name": "BUA",
                     "citytownclassification": "Classification",
                     "latitude": "lat",
                     "longitude": "lon",
                     "region_name": "Region",
                     "la_name": "Local Authority",
                     "lsoa_code": "LSOA code",
                     "constituency_name": "Constituency",
                     "population": "population",
                     "outputarea_code": "OA code"
                 }
                 var k = Object.keys(this_feature).filter(i=>(i != "type" && i != "geometry" && i != "Town" && i != "is_in" && i != "z_order" && i != "railway" && i != "pgroup" && i != "urban"))
                     for (var i = 0; i < k.length; i++) {
                         if (k[i] in this_feature) {
                             popupContent += '<br>' + lookup[k[i]] + ': ' + this_feature[k[i]];
                         }
                     }
                 layer.bindPopup(popupContent);
             }
             L.geoJSON(d, {
                 style: function(feature) {                     
                     switch (feature.geometry.type) {
                         case 'Point':
                             return {
                                 color: "red",
                                 radius: radius,
                                 weight: weight
                             };
                         default:
                             return {
                                 weight: weight
                             };
                     }
                 },
                 onEachFeature: onEachFeature,
                 pointToLayer: function(feature, latlng) {
                     return L.circleMarker(latlng, {
                         opacity: 1,
                         fillOpacity: 0.8
                     });
                 }
             }).addTo(map);
         });
        </script>
    </body>
 </html>
diff --git a/LICENSE b/LICENSE
 MIT License

 Copyright (c) 2020 Will Deakin

 Permission is hereby granted, free of charge, to any person obtaining a copy
 of this software and associated documentation files (the "Software"), to deal
 in the Software without restriction, including without limitation the rights
 to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
 copies of the Software, and to permit persons to whom the Software is
 furnished to do so, subject to the following conditions:

 The above copyright notice and this permission notice shall be included in all
 copies or substantial portions of the Software.

 THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
 IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
 FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
 AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
 LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
 OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
 SOFTWARE.

diff --git a/prepublish b/prepublish
 #!/bin/sh

 export PYTHONUNBUFFERED=1

 for i in 'download' 'work' 'shp'
 do
    if [ ! -d ${i} ]; then
        mkdir -p ${i}
    fi
 done

 URL=https://www.nrscotland.gov.uk/files/geography
 #FILENAME=output-area-2011-eor.zip
 FILENAME=output-area-2011-mhw.zip
 if [ ! -s download/${FILENAME} ]; then
    curl -L -o download/${FILENAME} ${URL}/${FILENAME}
 fi

 if [ ! -s work/OutputArea2011_MHW.shp ]; then
    (cd work; unzip ../download/${FILENAME})
 fi

 URL=https://opendata.arcgis.com/datasets
 #FILENAME=d74074ae6dec4de59fdcd2744fefc1f9_0.zip
 FILENAME=f79fc19485704ce68523d8d70d84a913_0.zip
 if [ ! -s download/${FILENAME} ]; then
    curl -L -o download/${FILENAME} ${URL}/${FILENAME}
 fi

 if [ ! -s work/Output_Areas__December_2011__Boundaries_EW_BGC.shp ]; then    
    (cd work; unzip ../download/${FILENAME})
 fi

 URL=https://researchbriefings.files.parliament.uk/documents/CBP-8322
 FILENAME=oa-classification-csv.csv
 if [ ! -s ${FILENAME} ]; then
    curl -L -o ${FILENAME} ${URL}/${FILENAME}
 fi

 ./town-boundaries.py
diff --git a/requirements.txt b/requirements.txt
 geopandas >= 0.9.0
 xlrd >= 2.0.1
diff --git a/town-boundaries.py b/town-boundaries.py
 #!/usr/bin/env python3

 import os
 from functools import partial
 from fiona.transform import transform_geom

 import pandas as pd
 import geopandas as gp
 from shapely.ops import transform

 # EPSG:4326 WG 84
 # EPSG:32630
 # EPSG:27700 OS GB36

 URBANTYPES = {'Large Town',
         'Large Town in Conurbation',
         'Core City (outside London)',
         'Village or small community in Conurbation',
         'Other City',
         'Small Town in Conurbation',
         'Small Town',
         'Medium Town',
         'Medium Town in Conurbation',
         'Core City (London)',
         'Village or Small Community in Conurbation'}

 def _set_precision(precision=6):
    def _precision(x, y, z=None):
        return tuple([round(i, precision) for i in [x, y, z] if i])
    return partial(transform, _precision)

 pd.set_option('display.max_columns', None)

 print('Load Output Area Data')
 TOWNDATA = pd.read_csv('oa-classification-csv.csv')

 TOWNDATA['name'] = TOWNDATA['bua_name']
 IDX1 = TOWNDATA['bua_name'] == 'None'
 TOWNDATA.loc[IDX1, 'name'] = TOWNDATA.loc[IDX1, 'la_name']
 for k in [' BUA in Conurbation', ' BUASD', ' BUA']:
    TOWNDATA['name'] = TOWNDATA['name'].str.replace(k, '')

 TOWNS = TOWNDATA.groupby(['bua_code', 'name', 'region_name']).filter(lambda v: v['population'].sum() > 1)
 TOWNS['Town'] = TOWNDATA['name']
 TOWNDATA = None

 IDX2 = TOWNS['citytownclassification'].isin(URBANTYPES)
 TOWNS['urban'] = 'N'
 TOWNS.loc[IDX2, 'urban'] = 'Y'

 print('Loaded Output Area Data')
 print('Load Scotland')
 SCDATA = gp.read_file('work/OutputArea2011_MHW.shp')
 CRS = SCDATA.crs.srs

 print('Aggregate Scotland')
 SC = SCDATA[['DataZone', 'geometry']].dissolve(by='DataZone', aggfunc='sum')

 SDX1 = TOWNS['lsoa_code'].str[:1] == 'S'
 G1 = TOWNS.loc[SDX1, ['lsoa_code']].join(SC['geometry'], on='lsoa_code')
 SC = None
 SCDATA = None
 print('Loaded Scotland')

 print('Load England and Wales')
 EW = gp.read_file('work/Output_Areas__December_2011__Boundaries_EW_BGC.shp')
 EW = EW.to_crs(CRS)
 EW = EW.set_index('OA11CD')
 G2 = TOWNS.loc[~SDX1, ['outputarea_code']].join(EW['geometry'], on='outputarea_code')
 EW = None
 print('Loaded England and Wales')

 print('Create GeoDataFrame')
 GEOMETRY = pd.concat([G2['geometry'], G1['geometry']])
 G1 = None
 G2 = None

 GF1 = gp.GeoDataFrame(TOWNS, geometry=GEOMETRY, crs=CRS)
 _precision = _set_precision(1)
 GF1['geometry'] = GF1['geometry'].apply(_precision)

 print('Write GeoDataFrame')
 GF1.to_file('all_boundaries.geojson', crs=CRS, driver='GeoJSON')
 GF1.to_file('shp/all_boundaries.shp', crs=CRS)

 GEOMETRY = None
 print('Aggregate by region, town, Built Up Area (BUA) and Middle Super-Output-Area')
 KEYS = ['region_name', 'Town', 'bua_code', 'msoa_code']
 GF2 = GF1[KEYS + ['geometry', 'population']].dissolve(by=KEYS, aggfunc=sum)
 #GF2 = GF2.to_crs('EPSG:4326')

 print('Aggregated data')
 FULLKEYS = KEYS + ['outputarea_code', 'lsoa_code', 'la_name', 'bua_name', 'constituency_name', 'citytownclassification', 'urban']
 DF1 = TOWNS[FULLKEYS].drop_duplicates(subset=KEYS).set_index(KEYS)

 _precision = _set_precision(5)
 POINTS = gp.GeoDataFrame(geometry=GF2.centroid, crs=CRS).to_crs('EPSG:4326')
 POINTS['geometry'] = POINTS['geometry'].apply(_precision)
 POINTS['data'] = POINTS['geometry'].apply(lambda v: [v.x, v.y])
 POINTS['longitude'], POINTS['latitude'] = zip(*POINTS.pop('data'))
 DF2 = POINTS[['longitude', 'latitude']]

 BOUNDARIES = GF2.join(DF1).join(DF2).reset_index().fillna('-')
 BOUNDARIES = BOUNDARIES.to_crs('EPSG:4326')
 for k in ['bua_code', 'bua_name', 'citytownclassification']:
    BOUNDARIES[k] = BOUNDARIES[k].str.replace('None', '-')

 #from matplotlib import pyplot as plt
 #fig = plt.figure()
 #ax1 = fig.add_subplot(111)

 print('Segment boundaries')
 INTERVAL = {'rural': 0, '5k': 5E3, '10k': 10E3, '50k': 50E3, '100k': 100E3, 'MC1': 1E9}
 BINS = pd.IntervalIndex.from_breaks(list(INTERVAL.values()), closed='left')
 KEYS = ['region_name', 'Town', 'bua_code']
 POPULATION = BOUNDARIES[KEYS + ['population']].groupby(KEYS).sum()
 CUT = pd.cut(POPULATION['population'].to_numpy(), BINS)
 CUT = CUT.rename_categories({i: j for i, j in zip(CUT.categories, INTERVAL.keys())})
 POPULATION['pgroup'] = CUT.tolist()
 BOUNDARIES = BOUNDARIES.set_index(KEYS).join(POPULATION['pgroup']).reset_index()
 IDX3 = BOUNDARIES['urban'] == 'N'
 BOUNDARIES.loc[IDX3, 'pgroup'] = 'rural'

 print('Write full boundaries')
 _precision = _set_precision(5)
 BOUNDARIES['geometry'] = BOUNDARIES['geometry'].apply(_precision)
 BOUNDARIES.to_file('shp/full_boundaries.shp', crs='epsg:4326')
 BOUNDARIES.to_file('full_boundaries.geojson', crs='epsg:4326', driver='GeoJSON')

 print('Write segmented boundaries')
 for k, g in BOUNDARIES.groupby('pgroup'):
    stub = 'town_boundaries_{}'.format(k)
    g.to_file('shp/{}.shp'.format(stub), crs='epsg:4326')
    g.to_file('{}.geojson'.format(stub), driver='GeoJSON', crs='epsg:4326')

 BOUNDARIES = BOUNDARIES.to_crs('epsg:27700')
 for k, g in BOUNDARIES.groupby('pgroup'):
    stub = 'town_boundaries_gb_{}'.format(k)
    g.to_file('shp/{}.shp'.format(stub), crs='epsg:27700')

diff --git a/town_boundaries.geojson b/town_boundaries.geojson
diff --git a/town_boundaries_100k.geojson b/town_boundaries_100k.geojson
diff --git a/town_boundaries_10k.geojson b/town_boundaries_10k.geojson
diff --git a/town_boundaries_50k.geojson b/town_boundaries_50k.geojson
diff --git a/town_boundaries_5k.geojson b/town_boundaries_5k.geojson
diff --git a/town_boundaries_rural.geojson b/town_boundaries_rural.geojson
	Released license:MIT
	height:780
	border:no

	Towns: Location of urbans areas > 10k population in Great Britain

	Digital boundaries and reference map data:
	Source: Office for National Statistics licensed under the Open Government Licence v.3.0.
	Source: National Records Scotland data © Crown copyright and database right 2020.
	Contains OS data © Crown copyright and database right 2020.
	<!DOCTYPE html>
	<html>
	<head>
	<title>Towns: Location in GB</title>
	<meta charset="utf-8"/>
	<meta name="viewport" content="width=device-width, initial-scale=1.0">
	<meta http-equiv="Cache-control" content="No-Cache">
	<link rel="shortcut icon" type="image/x-icon" href="docs/images/favicon.ico"/>
	<link rel="stylesheet" href="https://unpkg.com/[email protected]/dist/leaflet.css"
	integrity="sha512-xwE/Az9zrjBIphAcBb3F6JVqxf46+CDLwfLMHloNu6KEQCAWi6HcDUbeOfBIptF7tcCzusKFjFw2yuvEpDL9wQ==" crossorigin=""/>
	<script src="https://unpkg.com/[email protected]/dist/leaflet.js"
	integrity="sha512-gZwIG9x3wUXg2hdXF6+rVkLF/0Vi9U8D2Ntg4Ga5I5BZpVkVxlJWbSQtXPSiUTtC0TjtGOmxa1AJPuV0CPthew=="
	crossorigin=""></script>
	<script src="http://d3js.org/d3.v5.min.js"></script>
	<div id="map" style="width: 1340px; height: 780px"></div>
	</head>
	<body>
	<div id='map'></div>
	<script type="text/javascript">
	var radius = 3;
	var weight = 1;
	var linewidth = 2;
	var log2 = Math.log(2.0);
	var minZoom = 3;
	var maxZoom = 18;
	var map = L.map('map').setView([54.533, -2.53], 6);

	L.tileLayer("http://{s}.tile.openstreetmap.org/{z}/{x}/{y}.png", {
	maxZoom: maxZoom,
	minZoom: minZoom,
	attribution: 'Map data © <a href="https://www.openstreetmap.org/">OpenStreetMap</a> contributors, ' + '<a href="https://creativecommons.org/licenses/by-sa/2.0/">CC-BY-SA</a>, ' + 'Imagery © <a href="https://www.mapbox.com/">Mapbox</a>',
	id: 'osm.standard'
	}).addTo(map);

	d3.json("town_boundaries.geojson").then(function(d) {
	function onEachFeature(feature, layer) {
	var this_feature = feature.properties;
	var popupContent;
	if (this_feature.Town)
	popupContent = this_feature.Town;
	else
	popupContent = this_feature.railway;
	var lookup = {
	"bua_code": "BUA code",
	"bua_name": "BUA",
	"citytownclassification": "Classification",
	"latitude": "lat",
	"longitude": "lon",
	"region_name": "Region",
	"la_name": "Local Authority",
	"lsoa_code": "LSOA code",
	"constituency_name": "Constituency",
	"population": "population",
	"outputarea_code": "OA code"
	}
	var k = Object.keys(this_feature).filter(i=>(i != "type" && i != "geometry" && i != "Town" && i != "is_in" && i != "z_order" && i != "railway" && i != "pgroup" && i != "urban"))
	for (var i = 0; i < k.length; i++) {
	if (k[i] in this_feature) {
	popupContent += '<br>' + lookup[k[i]] + ': ' + this_feature[k[i]];
	}
	}
	layer.bindPopup(popupContent);
	}
	L.geoJSON(d, {
	style: function(feature) {
	switch (feature.geometry.type) {
	case 'Point':
	return {
	color: "red",
	radius: radius,
	weight: weight
	};
	default:
	return {
	weight: weight
	};
	}
	},
	onEachFeature: onEachFeature,
	pointToLayer: function(feature, latlng) {
	return L.circleMarker(latlng, {
	opacity: 1,
	fillOpacity: 0.8
	});
	}
	}).addTo(map);
	});
	</script>
	</body>
	</html>
	MIT License

	Copyright (c) 2020 Will Deakin

	Permission is hereby granted, free of charge, to any person obtaining a copy
	of this software and associated documentation files (the "Software"), to deal
	in the Software without restriction, including without limitation the rights
	to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
	copies of the Software, and to permit persons to whom the Software is
	furnished to do so, subject to the following conditions:

	The above copyright notice and this permission notice shall be included in all
	copies or substantial portions of the Software.

	THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
	IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
	FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
	AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
	LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
	OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
	SOFTWARE.
	#!/bin/sh

	export PYTHONUNBUFFERED=1

	for i in 'download' 'work' 'shp'
	do
	if [ ! -d ${i} ]; then
	mkdir -p ${i}
	fi
	done

	URL=https://www.nrscotland.gov.uk/files/geography
	#FILENAME=output-area-2011-eor.zip
	FILENAME=output-area-2011-mhw.zip
	if [ ! -s download/${FILENAME} ]; then
	curl -L -o download/${FILENAME} ${URL}/${FILENAME}
	fi

	if [ ! -s work/OutputArea2011_MHW.shp ]; then
	(cd work; unzip ../download/${FILENAME})
	fi

	URL=https://opendata.arcgis.com/datasets
	#FILENAME=d74074ae6dec4de59fdcd2744fefc1f9_0.zip
	FILENAME=f79fc19485704ce68523d8d70d84a913_0.zip
	if [ ! -s download/${FILENAME} ]; then
	curl -L -o download/${FILENAME} ${URL}/${FILENAME}
	fi

	if [ ! -s work/Output_Areas__December_2011__Boundaries_EW_BGC.shp ]; then
	(cd work; unzip ../download/${FILENAME})
	fi

	URL=https://researchbriefings.files.parliament.uk/documents/CBP-8322
	FILENAME=oa-classification-csv.csv
	if [ ! -s ${FILENAME} ]; then
	curl -L -o ${FILENAME} ${URL}/${FILENAME}
	fi

	./town-boundaries.py
	#!/usr/bin/env python3

	import os
	from functools import partial
	from fiona.transform import transform_geom

	import pandas as pd
	import geopandas as gp
	from shapely.ops import transform

	# EPSG:4326 WG 84
	# EPSG:32630
	# EPSG:27700 OS GB36

	URBANTYPES = {'Large Town',
	'Large Town in Conurbation',
	'Core City (outside London)',
	'Village or small community in Conurbation',
	'Other City',
	'Small Town in Conurbation',
	'Small Town',
	'Medium Town',
	'Medium Town in Conurbation',
	'Core City (London)',
	'Village or Small Community in Conurbation'}

	def _set_precision(precision=6):
	def _precision(x, y, z=None):
	return tuple([round(i, precision) for i in [x, y, z] if i])
	return partial(transform, _precision)

	pd.set_option('display.max_columns', None)

	print('Load Output Area Data')
	TOWNDATA = pd.read_csv('oa-classification-csv.csv')

	TOWNDATA['name'] = TOWNDATA['bua_name']
	IDX1 = TOWNDATA['bua_name'] == 'None'
	TOWNDATA.loc[IDX1, 'name'] = TOWNDATA.loc[IDX1, 'la_name']
	for k in [' BUA in Conurbation', ' BUASD', ' BUA']:
	TOWNDATA['name'] = TOWNDATA['name'].str.replace(k, '')

	TOWNS = TOWNDATA.groupby(['bua_code', 'name', 'region_name']).filter(lambda v: v['population'].sum() > 1)
	TOWNS['Town'] = TOWNDATA['name']
	TOWNDATA = None

	IDX2 = TOWNS['citytownclassification'].isin(URBANTYPES)
	TOWNS['urban'] = 'N'
	TOWNS.loc[IDX2, 'urban'] = 'Y'

	print('Loaded Output Area Data')
	print('Load Scotland')
	SCDATA = gp.read_file('work/OutputArea2011_MHW.shp')
	CRS = SCDATA.crs.srs

	print('Aggregate Scotland')
	SC = SCDATA[['DataZone', 'geometry']].dissolve(by='DataZone', aggfunc='sum')

	SDX1 = TOWNS['lsoa_code'].str[:1] == 'S'
	G1 = TOWNS.loc[SDX1, ['lsoa_code']].join(SC['geometry'], on='lsoa_code')
	SC = None
	SCDATA = None
	print('Loaded Scotland')

	print('Load England and Wales')
	EW = gp.read_file('work/Output_Areas__December_2011__Boundaries_EW_BGC.shp')
	EW = EW.to_crs(CRS)
	EW = EW.set_index('OA11CD')
	G2 = TOWNS.loc[~SDX1, ['outputarea_code']].join(EW['geometry'], on='outputarea_code')
	EW = None
	print('Loaded England and Wales')

	print('Create GeoDataFrame')
	GEOMETRY = pd.concat([G2['geometry'], G1['geometry']])
	G1 = None
	G2 = None

	GF1 = gp.GeoDataFrame(TOWNS, geometry=GEOMETRY, crs=CRS)
	_precision = _set_precision(1)
	GF1['geometry'] = GF1['geometry'].apply(_precision)

	print('Write GeoDataFrame')
	GF1.to_file('all_boundaries.geojson', crs=CRS, driver='GeoJSON')
	GF1.to_file('shp/all_boundaries.shp', crs=CRS)

	GEOMETRY = None
	print('Aggregate by region, town, Built Up Area (BUA) and Middle Super-Output-Area')
	KEYS = ['region_name', 'Town', 'bua_code', 'msoa_code']
	GF2 = GF1[KEYS + ['geometry', 'population']].dissolve(by=KEYS, aggfunc=sum)
	#GF2 = GF2.to_crs('EPSG:4326')

	print('Aggregated data')
	FULLKEYS = KEYS + ['outputarea_code', 'lsoa_code', 'la_name', 'bua_name', 'constituency_name', 'citytownclassification', 'urban']
	DF1 = TOWNS[FULLKEYS].drop_duplicates(subset=KEYS).set_index(KEYS)

	_precision = _set_precision(5)
	POINTS = gp.GeoDataFrame(geometry=GF2.centroid, crs=CRS).to_crs('EPSG:4326')
	POINTS['geometry'] = POINTS['geometry'].apply(_precision)
	POINTS['data'] = POINTS['geometry'].apply(lambda v: [v.x, v.y])
	POINTS['longitude'], POINTS['latitude'] = zip(*POINTS.pop('data'))
	DF2 = POINTS[['longitude', 'latitude']]

	BOUNDARIES = GF2.join(DF1).join(DF2).reset_index().fillna('-')
	BOUNDARIES = BOUNDARIES.to_crs('EPSG:4326')
	for k in ['bua_code', 'bua_name', 'citytownclassification']:
	BOUNDARIES[k] = BOUNDARIES[k].str.replace('None', '-')

	#from matplotlib import pyplot as plt
	#fig = plt.figure()
	#ax1 = fig.add_subplot(111)

	print('Segment boundaries')
	INTERVAL = {'rural': 0, '5k': 5E3, '10k': 10E3, '50k': 50E3, '100k': 100E3, 'MC1': 1E9}
	BINS = pd.IntervalIndex.from_breaks(list(INTERVAL.values()), closed='left')
	KEYS = ['region_name', 'Town', 'bua_code']
	POPULATION = BOUNDARIES[KEYS + ['population']].groupby(KEYS).sum()
	CUT = pd.cut(POPULATION['population'].to_numpy(), BINS)
	CUT = CUT.rename_categories({i: j for i, j in zip(CUT.categories, INTERVAL.keys())})
	POPULATION['pgroup'] = CUT.tolist()
	BOUNDARIES = BOUNDARIES.set_index(KEYS).join(POPULATION['pgroup']).reset_index()
	IDX3 = BOUNDARIES['urban'] == 'N'
	BOUNDARIES.loc[IDX3, 'pgroup'] = 'rural'

	print('Write full boundaries')
	_precision = _set_precision(5)
	BOUNDARIES['geometry'] = BOUNDARIES['geometry'].apply(_precision)
	BOUNDARIES.to_file('shp/full_boundaries.shp', crs='epsg:4326')
	BOUNDARIES.to_file('full_boundaries.geojson', crs='epsg:4326', driver='GeoJSON')

	print('Write segmented boundaries')
	for k, g in BOUNDARIES.groupby('pgroup'):
	stub = 'town_boundaries_{}'.format(k)
	g.to_file('shp/{}.shp'.format(stub), crs='epsg:4326')
	g.to_file('{}.geojson'.format(stub), driver='GeoJSON', crs='epsg:4326')

	BOUNDARIES = BOUNDARIES.to_crs('epsg:27700')
	for k, g in BOUNDARIES.groupby('pgroup'):
	stub = 'town_boundaries_gb_{}'.format(k)
	g.to_file('shp/{}.shp'.format(stub), crs='epsg:27700')