[PYTHON] Try to display the railway data of national land numerical information in 3D
Let's display the railway lines that can be obtained from the railway data of the national land numerical information in 3D.
demo http://needtec.sakura.ne.jp/threemap/railroad3d_example.html
** Client-side source ** https://github.com/mima3/threemap
** Server-side source ** https://github.com/mima3/kokudo
Server-side processing
The server-side process returns the stored railroad data of the national land numerical information under the requested conditions. At this time, altitude data is added.
Import railroad data into SpatiaLite from national land numerical information
Import Railway Information of the national land numerical information into SpatiaLite.
The operation of spatialite is as follows. https://github.com/mima3/kokudo/blob/master/kokudo_db.py
To import the shapefile of the railway data of national land numerical information into the database, execute the following command.
python import_railroad_section.py C:\tool\spatialite\mod_spatialite-4.2.0-win-x86\mod_spatialite.dll test.sqlite original_data\N02-13\N02-13_RailroadSection.shp
Allow GeoJson to be retrieved via HTTP.
Allows you to get Geojson for a specified route from a database created via HTTP using Bottle.
** Call example ** http://needtec.sakura.ne.jp/kokudo/json/get_railroad_section?operationCompany=%E6%9D%B1%E4%BA%AC%E6%80%A5%E8%A1%8C%E9%9B%BB%E9%89%84
In this example, all lines of Tokyu Corporation are acquired.
Find the elevation using the Geographical Survey Map elevation API
The railway data of the national land numerical information has longitude and latitude, but not altitude. There is no point in displaying it in 3D.
Therefore, we will make it possible to obtain the altitude from the longitude and latitude. To do this, use the Geographical Survey Map elevation API.
http://portal.cyberjapan.jp/help/development/api.html
Example of use: http://cyberjapandata2.gsi.go.jp/general/dem/scripts/getelevation.php?lon=140.08531&lat=36.103543&outtype=JSON
result
{"elevation":25.3,"hsrc":"5m\uff08\u30ec\u30fc\u30b6\uff09"}
However, it is useless to execute this API every time, so make sure to cache the contents once read in the database.
py:https://github.com/mima3/kokudo/blob/master/gsi_api.py
# -*- coding: utf-8 -*-
import urllib
import urllib2
from peewee import *
from playhouse.sqlite_ext import SqliteExtDatabase
import json
database_proxy = Proxy() # Create a proxy for our db.
def get_elevation_by_api(long, lat):
"""
Obtaining the altitude value
http://portal.cyberjapan.jp/help/development/api.html
"""
url = ('http://cyberjapandata2.gsi.go.jp/general/dem/scripts/getelevation.php?lon=%f&lat=%f&outtype=JSON' % (long, lat))
req = urllib2.Request(url)
opener = urllib2.build_opener()
conn = opener.open(req)
cont = conn.read()
ret = json.loads(cont)
return ret
def str_isfloat(str):
try:
float(str)
return True
except ValueError:
return False
class ElevationCache(Model):
"""
Elevation cache table
"""
lat = FloatField(index=True)
long = FloatField(index=True)
hsrc = TextField(index=True)
elevation = FloatField(null=True)
class Meta:
database = database_proxy
def connect(path):
db = SqliteExtDatabase(path)
database_proxy.initialize(db)
def setup(path):
connect(path)
database_proxy.create_tables([ElevationCache], True)
def get_elevation(long, lat):
try:
ret = ElevationCache.get((ElevationCache.long==long) & (ElevationCache.lat==lat))
return {'elevation': ret.elevation, 'hsrc': ret.hsrc}
except ElevationCache.DoesNotExist:
ret = get_elevation_by_api(long, lat)
elevation = ret['elevation']
if not str_isfloat(elevation):
elevation = None
ElevationCache.create(
long = long,
lat = lat,
elevation = elevation,
hsrc = ret['hsrc']
)
return ret
class GsiConvertError(Exception):
def __init__(self, value):
self.value = value
def __str__(self):
return repr(self.value)
def convert_geojson(json):
for feature in json['features']:
if feature['geometry']['type'] == 'LineString':
prop = {}
start = feature['geometry']['coordinates'][0]
end = feature['geometry']['coordinates'][len(feature['geometry']['coordinates'])-1]
start_elevation = get_elevation(start[0], start[1])
end_elevation = get_elevation(end[0], end[1])
feature['properties']['start_elevation'] = start_elevation['elevation']
feature['properties']['end_elevation'] = end_elevation['elevation']
else:
raise GsiConvertError('unexpected feature type')
return json
if __name__ == '__main__':
setup('elevation_cache.sqlite')
#print get_elevation_by_api(133, 39)
#print get_elevation_by_api(139.766084, 35.681382)
print get_elevation(133, 39)
print get_elevation(139.766084, 35.681382)
with open('get_railroad_section.geojson' , 'rb') as f:
cont = f.read()
print convert_geojson(json.loads(cont))
The get_elevation function finds the altitude from longitude and latitude. At this time, if the value is already stored in the DB, it is returned, if it is not stored, the elevation API is executed, and the result is stored in the DB and returned.
The convert_geojson function assigns elevation to GeoJson for LineString. Gets the elevation for the start and end points of the line and stores the results in the properties as start_elevation, end_elevation.
Get GeoJSON of railway data with elevation data.
Example of use http://needtec.sakura.ne.jp/kokudo/json/get_railroad_section?operationCompany=%E6%9D%B1%E4%BA%AC%E6%80%A5%E8%A1%8C%E9%9B%BB%E9%89%84&embed_elevation=True
Elevation can be obtained by adding embed_elevation.
** Acquisition result **
{
"type": "FeatureCollection",
"features": [
{
"geometry": {
"type": "LineString",
"coordinates": [[139.48677, 35.55760999999999], [139.4865599999999, 35.55839]]
},
"type": "Feature",
"properties": {
"operationCompany": "\u6771\u4eac\u6025\u884c\u96fb\u9244",
"end_elevation": 36.7,
"serviceProviderType": "4",
"railwayLineName": "\u3053\u3069\u3082\u306e\u56fd\u7dda",
"railwayType": "12",
"start_elevation": 35.4
}
},//Abbreviation
]
}
Client-side processing
The processing on the client side is as follows.
-
- Give the acquired GeoJSON z-axis. At this time, the start point and the end point are the altitudes specified in properties, and the other points are the ratio considering the start point and the end point.
- Lines that can be combined are combined
-
- Drawing routes using TubeGeometry in three.js
Give the acquired GeoJSON z-axis
Add elevation to the z-axis for all coordinates using start_elevation and end_elevation in properties as shown below.
function expendElevation(features) {
for (var i = 0; i < features.length; ++i) {
var feature = features[i];
var end_elevation = feature.properties.end_elevation;
var start_elevation = feature.properties.start_elevation;
var per_elevation = (end_elevation - start_elevation) / feature.geometry.coordinates.length;
for (var j = 0; j < feature.geometry.coordinates.length; ++j) {
feature.geometry.coordinates[j].push(start_elevation + (j * per_elevation));
}
}
return features;
}
Lines that can be combined are combined
For example, if the start point of one feature and the end point of another feature are equal, they are combined and regarded as one line.
function compressionLine(features) {
var before = features.length;
for (var i = 0; i < features.length; ++i) {
for (var j = features.length -1; i < j; --j) {
var f1Start = features[i].geometry.coordinates[0];
var f1End = features[i].geometry.coordinates[features[i].geometry.coordinates.length-1];
var f2Start = features[j].geometry.coordinates[0];
var f2End = features[j].geometry.coordinates[features[j].geometry.coordinates.length-1];
//If the start point of f1 coincides with the end point of f2, then f2 precedes f1.
if (f1Start[0] == f2End[0] && f1Start[1] == f2End[1]) {
features[i].geometry.coordinates = features[j].geometry.coordinates.concat(features[i].geometry.coordinates);
features.splice(j, 1);
break;
}
//If the end point of f1 coincides with the start point of f2, then there is f2 after f1
if (f1End[0] == f2Start[0] && f1End[1] == f2Start[1]) {
features[i].geometry.coordinates = features[i].geometry.coordinates.concat(features[j].geometry.coordinates);
features.splice(j, 1);
break;
}
}
}
if (features.length == before) {
return features;
}
return compressionLine(features);
}
Draw routes using Tube Geometry
Draw a route using TubeGeometry. Determine the position of mesh.position with the coordinates of the start point as the offset. For other points, when adding them as TubeGeometry paths, give them relative coordinates to the starting point.
function createRailroad(geodata) {
console.log(geodata.length);
geodata = expendElevation(geodata);
geodata = compressionLine(geodata);
console.log(geodata.length);
var scaleElevation = 100;
for (var i = 0 ; i < geodata.length ; i++) {
var lineList = [];
var geoFeature = geodata[i];
var baseX;
var baseY;
for (var j = 0; j < geoFeature.geometry.coordinates.length; ++j) {
var pt = reverseProjection([
geoFeature.geometry.coordinates[j][0],
geoFeature.geometry.coordinates[j][1]
]);
if (j ==0) {
baseX = pt[0];
baseY = pt[1];
lineList.push(new THREE.Vector3(0, 0, geoFeature.geometry.coordinates[j][2] / scaleElevation));
} else {
lineList.push(new THREE.Vector3(pt[0] - baseX, pt[1] - baseY, geoFeature.geometry.coordinates[j][2] /scaleElevation));
}
}
var spline = new THREE.SplineCurve3(lineList);
var tubeGeo = new THREE.TubeGeometry(spline, 32, 0.03, 8, false);
var mesh = new THREE.Mesh(
tubeGeo,
new THREE.MeshLambertMaterial( {
color: 0xff0000,
transparent: true,
opacity: 0.9
})
);
mesh.position.set(baseX, baseY, 0.1);
scene.add(mesh);
}
}
At this time, the reverseProjection used to get the coordinates from the longitude and latitude uses the projection of d3.geo.path () as follows.
Contents of reverse Projection
//Since it is upside down, it is stuffy
var reverseProjection = function(x, y) {
pt = projection(x, y);
pt[1] *= -1;
return pt;
};
//Create a function to convert geoJSON data to a path
var path = d3.geo.path().projection(reverseProjection);
Summary
Geojson can be created dynamically by storing the railway data of national land numerical information in spatialite.
If you use the Geographical Survey Institute's elevation API, you can get the elevation from latitude and longitude. At this time, it is better to cache the result.
If you use three.js, you can easily perform 3D representation, and at this time, if you use d3 projection, you can easily convert longitude and latitude.
With these, you can express railway lines in 3D. ~~ But the subway and the Shinkansen can't handle this program, it's not good, you can't use it! ~~
Recommended Posts