Skip to content

Instantly share code, notes, and snippets.

@anisotropi4
Last active June 21, 2024 11:40
Show Gist options
  • Save anisotropi4/7ba18de4bfe15d9fea2bb5139dbc30e1 to your computer and use it in GitHub Desktop.
Save anisotropi4/7ba18de4bfe15d9fea2bb5139dbc30e1 to your computer and use it in GitHub Desktop.
OSM Rail Network in Great Britain
Released license:MIT
height:780
border:no
Trunk Roads: OSM locations of trunk roads in the GB
OpenStreetMap® open data is © OpenStreetMap contributors and used under the Open Street Map license https://www.openstreetmap.org/copyright
Open Data Institute Leeds scripts and configuration are used © Creative Commons BY 4.0 ODI Leeds 2018 https://creativecommons.org/licenses/by/4.0/

OSM Rail Network in Great Britain

Using OSMnx (Python for Street Maps) to download Open Street Map (OSM) rail-data using the Overpass API

Creating the datafiles and associate GeoJSON format report

Once the dependencies are met run the script to gather the data:

$ ./prepublish

This will download and create ESRI Shape files in the shp directory for all Towns as well as all locations in Great Britian, as well as a Leaflet JavaScript visualisation using the Town data in a GeoJSON format.

Dependencies

These are environment and project dependencies.

python dependencies

For ease of use manage package python packages dependencies with a local virtual environment venv:

$ virtualenv venv
$ source venv/bin/activate
$ pip install osmnx

If missing install the python virtualenv package

$ sudo apt install virtualenv

Data and Map Tile License

Open Street Map data is licensed under the Open Data Commons Open Database License (ODbL)

Open Street Map tiles are licenced by CC BY-SA (c) OpenStreetMap contributors

{"type": "Polygon", "coordinates": [[[1.98396, 51.3377], [0.058076, 50.10375], [-1.246112, 49.89032], [-2.027284, 49.77612], [-2.0209, 49.76352], [-2.03199, 49.76171], [-1.808598, 49.10263], [-1.835368, 49.00443], [-3.940324, 49.15592], [-7.02478, 49.70097], [-5.441616, 52.78017], [-5.206178, 53.77268], [-5.488813, 54.8651], [-6.208707, 55.33747], [-6.604158, 55.43057], [-7.148041, 55.66455], [-14.86751, 57.46894], [-14.9907, 57.68017], [-12.08287, 58.53333], [-1.637515, 61.13564], [-0.2670263, 61.103], [0.1278458, 59.76591], [0.097816, 59.42198], [0.799897, 55.79959], [1.702395, 54.43427], [2.25, 52.58], [1.98396, 51.3377]]]}
<!DOCTYPE html>
<html>
<head>
<title>Fizzy Knitting: OSM locations of rail electrification in the UK</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/leaflet@1.6.0/dist/leaflet.css"
integrity="sha512-xwE/Az9zrjBIphAcBb3F6JVqxf46+CDLwfLMHloNu6KEQCAWi6HcDUbeOfBIptF7tcCzusKFjFw2yuvEpDL9wQ==" crossorigin=""/>
<script src="https://unpkg.com/leaflet@1.6.0/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);
//var voltages = ["25000", "1500", "750", "630", "600", "550"]
var voltages = ["660", "750", "1000", "1500", "3000", "6000", "11000", "15000", "20000", "25000", "50000"];
var colours = ["#4300BD","#5A00B2","#6F00A7","#86009C","#FFCC00","#FFAF00","#FF9400","#FF7D00","#FF6600","#FF4500","#F00000","#D00000"];
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("output-all.json").then(function(d) {
function onEachFeature(feature, layer) {
var this_feature = feature.properties;
var popupContent;
if (this_feature.name)
popupContent = this_feature.name;
else
popupContent = this_feature.railway;
var lookup = {
"osm_id": "OSM ID",
"place": "Place",
"natural": "Natural",
"place": "Place",
"highway": "Road",
"ref": "CRS",
"admin_level": "admin",
"boundary": "boundary",
"voltage": "voltage",
"frequency": "frequency",
"electrified": "electrified"
}
var k = Object.keys(this_feature).filter(i=>(i != "type" && i != "geometry" && i != "name" && i != "is_in" && i != "z_order" && i != "railway"))
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);
}
function get_colour(p) {
var colours = {
contact_line: {
550: '#e69c65',
600: '#bd5d16',
750: '#ffce36',
1500: '#ff893f',
25000: '#ff3818',
0: '#9ca3a9',
},
rail: {
550: '#b3e297',
630: '#7bc643',
750: '#09843b',
},
'4th_rail': {
630: '#0a3797',
750: '#11a7e0',
},
'contact_line;rail': {
'25000;750': '#fe85a8',
},
'contact_line;4th_rail': {
'25000;750': '#9d176s',
},
no: '#2c2e35',
fail: '#ff00ff',
}
var e = p.electrified;
var f = p.frequency;
var v = p.voltage;
var r = p.railway
if (v && Number(v) < 550) v = '550';
if (v && Number(v) > 25000) return colours['fail'];
if (e == 'yes') e = 'contact_line';
if (e == '3rd_rail') e = 'rail';
if (! e) e = 'no';
if (e in colours) return colours[e][v] || colours[e];
return colours['fail'];
}
L.geoJSON(d, {
style: function(feature) {
var c = get_colour(feature.properties);
switch (feature.geometry.type) {
case 'Point':
return {
color: c,
radius: radius,
weight: weight
};
case 'LineString':
return {
color: c,
weight: linewidth
};
case 'MultiLineString':
return {
color: c,
weight: linewidth
};
case 'MultiPolygon':
return {
weight: linewidth,
opacity: 1.0,
color: "Green",
fillColor: "Green",
fillOpacity: 0.2
};
case 'GeometryCollection':
return {
radius: radius,
weight: linewidth,
opacity: 1.0,
color: "Green",
fillColor: "GreenYellow",
fillOpacity: 0.2
};
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.
#!/usr/bin/env python
import json
import osmnx as ox
from osmnx.utils import log
from osmnx.downloader import _osm_network_download
from functools import partial
from itertools import groupby
import argparse as args
from shapely.geometry import shape, Point, LineString
from shapely.ops import transform
import geopandas as gp
ox.config(use_cache=True, log_console=True)
try:
import pandas as pd
pd.set_option('display.max_columns', None)
except ImportError:
pass
PARSER = args.ArgumentParser(description='Download OSM rail data from OverpassAPI ')
PARSER.add_argument('--shpfiles', dest='shp', action='store_true',
default=False,
help='create ESRI ShapeFile data (default False)')
ARGS = PARSER.parse_args()
SHP = ARGS.shp
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)
def get_node(e):
node = {'geometry': Point(e['lon'], e['lat'])}
node['osmid'] = e['id']
if 'tags' in e:
for t in ox.settings.useful_tags_node:
if t in e['tags']:
node['tags'] = True
node[t] = e['tags'][t]
return node
def get_path(e):
path = {}
path['osmid'] = e['id']
# remove any consecutive duplicate elements in the list of nodes
g_list = groupby(e['nodes'])
path['nodes'] = [g[0] for g in g_list]
if 'tags' in e:
for t in ox.settings.useful_tags_path:
if t in e['tags']:
path[t] = e['tags'][t]
return path
def get_linestring(v):
return LineString([NODES[i]['geometry'] for i in v])
with open('great-britain.json', 'r') as fin:
POLYGON = shape(json.load(fin))
ox.settings.useful_tags_node = ['lon', 'timestamp', 'user', 'lat', 'ref', 'changeset', 'id', 'version', 'uid', 'ref:tiploc', 'name', 'electrified', 'frequency', 'voltage', 'railway']
ox.settings.useful_tags_path = ['bridge', 'tunnel', 'width', 'id', 'maxspeed', 'junction', 'name', 'version', 'uid', 'changeset', 'landuse', 'timestamp', 'user', 'ref', 'oneway', 'ref:tiploc', 'name', 'electrified', 'frequency', 'voltage', 'railway']
DATA = _osm_network_download(POLYGON, 'way["railway"]', '["railway"!~"highway|cycleway|footway|path|pedestrian|steps|corridor|elevator|escalator|proposed|construction|bridleway|abandoned|platform|raceway|service"]["railway"~"rail|subway|light_rail|tram|narrow_gauge"]')
log('Create node and path data')
NODES = {}
PATHS = {}
for osm_data in DATA:
for e in osm_data['elements']:
if e['type'] == 'node':
key = e['id']
NODES[key] = get_node(e)
if e['type'] == 'way': #osm calls network paths 'ways'
key = e['id']
PATHS[key] = get_path(e)
log('Create GeoPandas dataframe')
GF1 = gp.GeoDataFrame.from_dict(PATHS, orient='index')
GF1['geometry'] = GF1['nodes'].apply(get_linestring)
GF1 = GF1.drop('nodes', axis=1)
GF1.crs = 'EPSG:4326'
GF1['type'] = 'way'
IDX1 = GF1.within(POLYGON) | GF1.intersects(POLYGON)
GF1['location'] = 'GB'
GF1.loc[~IDX1, 'location'] = 'IE'
if SHP:
log('Output shapefiles')
GF2 = gp.GeoDataFrame.from_dict(data=NODES, orient='index')
GF2.crs = 'EPSG:4326'
GF2 = GF2.loc[GF2['tags'].dropna().index].drop(columns='tags', axis=1)
GF2['type'] = 'node'
IDX2 = GF2.within(POLYGON) | GF2.intersects(POLYGON)
GF2['location'] = 'GB'
GF2.loc[~IDX2, 'location'] = 'IE'
GF1.fillna('').to_file('shp/full_ways.shp', crs='epsg:4326')
GF2.fillna('').to_file('shp/full_nodes.shp', crs='epsg:4326')
log('Trim dataframe')
GF1 = GF1[IDX1].drop('location', axis=1)
log('Create GeoJSON')
OUTPUT = GF1.__geo_interface__
log('Clean up GeoJSON')
for i in OUTPUT['features']:
i['properties'] = {k: v for k, v in i['properties'].items()
if v is not None}
i.pop('bbox', None)
log('Write output file')
with open('output-all.json', 'w') as fout:
json.dump(OUTPUT, fout)
Display the source blob
Display the rendered blob
Raw
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
View raw

(Sorry about that, but we can’t show files that are this big right now.)

@anisotropi4
Copy link
Author

You're right. This was a previous internal osmnx module API call. I've rewritten the code here using the updated module. As it's broken I'll updated it this weekend.

@kevin-blaugrana
Copy link

You're right. This was a previous internal osmnx module API call. I've rewritten the code here using the updated module. As it's broken I'll updated it this weekend.

Sorry to bother you, but in your data folder (in the magpie repo), isn't the 'britain-simple.geojson' file a simple Polygon?
I can't figure out how you get all those lines in that geojson...

Could you explain to me the steps to follow? I'd like to do a similar job for France.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment