Orthographic projection to SVG

gistfile1.py

import csv, re
from hashlib import sha1
from math import cos, radians, sin, atan, degrees, pi, sqrt, hypot

csv.field_size_limit(13107200)

paths = dict()

with open('tm-world-borders-03.csv', 'r') as csvfile:
	csvreader = csv.reader(csvfile, delimiter=',', quotechar='"')
	csvreader.next()
	for a in csvreader:
		paths[a[5]] = a[0]


def project(lat, lng):
	lat = radians(lat)
	lng = radians(lng)
	lat0 = radians(-37)
	lng0 = radians(174.75)
	x = cos(lat) * sin(lng - lng0)
	y = cos(lat0) * sin(lat) - sin(lat0) * cos(lat) * cos(lng - lng0)
	c = sin(lat0) * sin(lat) + cos(lat0) * cos(lat) * cos(lng - lng0)
	return (x, y, c < 0)


foreground = dict()
background = dict()
fgsweep = dict()
fgsweep2 = dict()
bgsweep = dict()
fgshadow = dict()
bgshadow = dict()
bgshadow2 = dict()
for (name, r) in paths.iteritems():
	foreground[name] = list()
	background[name] = list()
	fgsweep[name] = list()
	fgsweep2[name] = list()
	bgsweep[name] = list()
	fgshadow[name] = list()
	bgshadow[name] = list()
	bgshadow2[name] = list()
	for p in r.split(')),(('):
		path = list()
		previous = None
		crossed = False

		backgroundparts = list()
		foregroundparts = list()

		for pt in re.findall('((-?\d+\.\d+) (-?\d+\.\d+))+', p):
			(x, y, c) = project(float(pt[2]), float(pt[1]))
			if not c:
				angle = atan(y / x)
				distance = sqrt(x * x + y * y)
				if x < 0:
					if y > 0:
						angle += pi
					else:
						angle -= pi
				angle = int(round(degrees(angle)))

			if previous is not None and previous != c:
				if previous:
					backgroundparts.append(path)
				else:
					foregroundparts.append(path)
				path = list()
				crossed = True
			previous = c
			path.append('%.2f,%.2f' % (x * 1000 + 1000, y * -1000 + 1000))

		if previous:
			backgroundparts.append(path)
		else:
			foregroundparts.append(path)

		if len(backgroundparts) > 0:
			try:
				for part in backgroundparts:
					minx = 2000
					maxx = -1000
					minxpoint = None
					maxxpoint = None
					for point in part:
						(x, y) = point.split(',')
						x = float(x) - 1000
						y = 1000 - float(y)
						angle = atan(y / x)
						if angle < 0 and maxx > 1.3:
							angle += pi * 2
						if angle > 0 and minx < -1.3:
							angle -= pi * 2
						if angle < minx:
							minx = angle
							minxpoint = point
						if angle > maxx:
							maxx = angle
							maxxpoint = point

					(x, y) = minxpoint.split(',')
					x = float(x) - 1000
					y = 1000 - float(y)
					scale = 1000 / hypot(x, y)

					bgshadow[name].append('1000,1000 %s %s' % (minxpoint, maxxpoint))
					minxshadowx = x * scale + 1000
					minxshadowy = 1000 - y * scale

					(x, y) = maxxpoint.split(',')
					x = float(x) - 1000
					y = 1000 - float(y)
					scale = 1000 / hypot(x, y)
					maxxshadowx = x * scale + 1000
					maxxshadowy = 1000 - y * scale

					bgsweep[name].append('%s %s %s,%s A 1000,1000 0 0 1 %s,%s' % (
						minxpoint, maxxpoint,
						maxxshadowx, maxxshadowy,
						minxshadowx, minxshadowy
					))
					fgsweep2[name].append('1000,1000 %s,%s A 1000,1000 0 0 1 %s,%s' % (
						maxxshadowx, maxxshadowy,
						minxshadowx, minxshadowy
					))
			except:
				pass

			if not crossed:
				background[name].append(" L ".join(backgroundparts[0]))
			else:
				backgroundstr = " L ".join(backgroundparts[0])
				for part in backgroundparts[1:]:
					first = part[0]
					arc = " A 1000,1000 0 0 0 %s" % first
					backgroundstr = backgroundstr + arc + " L " + " L ".join(part)

				first = backgroundparts[0][0]
				arc = " A 1000,1000 0 0 0 %s" % first
				background[name].append(backgroundstr + arc)

		if len(foregroundparts) > 0:
			try:
				for part in foregroundparts:
					minx = 2000
					maxx = -1000
					minxpoint = None
					maxxpoint = None
					for point in part:
						(x, y) = point.split(',')
						x = float(x) - 1000
						y = 1000 - float(y)
						if x == 0:
							continue
						angle = atan(y / x)
						if angle < 0 and maxx > 1.3:
							angle += pi * 2
						if angle > 0 and minx < -1.3:
							angle -= pi * 2
						if angle < minx:
							minx = angle
							minxpoint = point
						if angle > maxx:
							maxx = angle
							maxxpoint = point

					(x, y) = minxpoint.split(',')
					x = float(x) - 1000
					y = 1000 - float(y)
					scale = 1100 / hypot(x, y)
					if scale > 10:
						continue

					fgsweep[name].append('1000,1000 %s %s' % (minxpoint, maxxpoint))
					minxshadowx = x * scale + 1000
					minxshadowy = 1000 - y * scale

					(x, y) = maxxpoint.split(',')
					x = float(x) - 1000
					y = 1000 - float(y)
					scale = 1100 / hypot(x, y)
					maxxshadowx = x * scale + 1000
					maxxshadowy = 1000 - y * scale

					fgshadow[name].append('%s %s %s,%s A 1100,1100 0 0 1 %s,%s' % (
						minxpoint, maxxpoint,
						maxxshadowx, maxxshadowy,
						minxshadowx, minxshadowy
					))
					bgshadow2[name].append('1000,1000 %s,%s A 1100,1100 0 0 1 %s,%s' % (
						maxxshadowx, maxxshadowy,
						minxshadowx, minxshadowy
					))
			except:
				pass

			if not crossed:
				foreground[name].append(" L ".join(foregroundparts[0]))
			else:
				foregroundstr = " L ".join(foregroundparts[0])
				for part in foregroundparts[1:]:
					first = part[0]
					arc = " A 1000,1000 0 0 1 %s" % first
					foregroundstr = foregroundstr + arc + " L " + " L ".join(part)

				first = foregroundparts[0][0]
				arc = " A 1000,1000 0 0 1 %s" % first
				foreground[name].append(foregroundstr + arc)
print '<?xml version="1.0" ?>'
print '''<svg width="4400" height="2200" viewBox="-100 -100 4400 2200" xmlns="http://www.w3.org/2000/svg"
	xmlns:inkscape="http://www.inkscape.org/namespaces/inkscape">'''

print '\t<g inkscape:label="bgsweep" inkscape:groupmode="layer" '\
	'transform="matrix(-1,0,0,1,4200,0)">'
for (name, paths) in bgsweep.iteritems():
	if len(paths) == 0:
		continue
	print '\t\t<g id="%s_bgsweep" fill="#%s">' % \
		(name, sha1(name).hexdigest()[0:6])
	for path in paths:
		print '\t\t\t<path d="M %s Z" />' % path
	print '\t\t</g>'
print '\t</g>'
print '\t<g inkscape:label="bgshadow" inkscape:groupmode="layer" '\
	'transform="matrix(-1,0,0,1,4200,0)">'
for (name, paths) in bgshadow.iteritems():
	if len(paths) == 0:
		continue
	print '\t\t<g id="%s_bgshadow" fill="#eee">' % name
	for path in paths:
		print '\t\t\t<path d="M %s Z" />' % path
	print '\t\t</g>'
print '\t</g>'
print '\t<g inkscape:label="bgshadow2" inkscape:groupmode="layer" '\
	'transform="matrix(-1,0,0,1,4200,0)">'
for (name, paths) in bgshadow2.iteritems():
	if len(paths) == 0:
		continue
	print '\t\t<g id="%s_bgshadow2" fill="#eee">' % name
	for path in paths:
		print '\t\t\t<path d="M %s Z" />' % path
	print '\t\t</g>'
print '\t</g>'
print '\t<g inkscape:label="background" inkscape:groupmode="layer" '\
	'transform="matrix(-1,0,0,1,4200,0)">'
for (name, paths) in background.iteritems():
	if len(paths) == 0:
		continue
	print '\t\t<g id="%s_bg">' % name
	for path in paths:
		print '\t\t\t<path d="M %s Z" fill="black" />' % path
	print '\t\t</g>'
print '\t</g>'
print '\t<g inkscape:label="fgsweep2" inkscape:groupmode="layer">'
for (name, paths) in fgsweep2.iteritems():
	if len(paths) == 0:
		continue
	print '\t\t<g id="%s_fgsweep2" fill="#%s">' % \
		(name, sha1(name).hexdigest()[0:6])
	for path in paths:
		print '\t\t\t<path d="M %s Z" />' % path
	print '\t\t</g>'
print '\t</g>'
print '\t<g inkscape:label="fgsweep" inkscape:groupmode="layer">'
for (name, paths) in fgsweep.iteritems():
	if len(paths) == 0:
		continue
	print '\t\t<g id="%s_fgsweep" fill="#%s">' % \
		(name, sha1(name).hexdigest()[0:6])
	for path in paths:
		print '\t\t\t<path d="M %s Z" />' % path
	print '\t\t</g>'
print '\t</g>'
print '\t<g inkscape:label="fgshadow" inkscape:groupmode="layer">'
for (name, paths) in fgshadow.iteritems():
	if len(paths) == 0:
		continue
	print '\t\t<g id="%s_fgshadow" fill="#eee">' % name
	for path in paths:
		print '\t\t\t<path d="M %s Z" />' % path
	print '\t\t</g>'
print '\t</g>'
print '\t<g inkscape:label="foreground" inkscape:groupmode="layer">'
for (name, paths) in foreground.iteritems():
	if len(paths) == 0:
		continue
	print '\t\t<g id="%s_fg" fill="red">' % name
	for path in paths:
		print '\t\t\t<path d="M %s Z" />' % path
	print '\t\t</g>'
print '\t</g>'
print '</svg>'