Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
Generator for Dorling cartograms
"""
Generator for packed circle cartograms
"""
import proj, gisutils
class Cartogram:
def loadCSV(self, url, key='id', value='val', lon='lon', lat='lat'):
import csv
doc = csv.reader(open(url))
head = None
circles = []
for row in doc:
if not head:
head = row
else:
circles.append(Circle(row[head.index(lon)], row[head.index(lat)], row[head.index(key)], row[head.index(value)]))
self.circles = circles
self.computeRadii()
def computeRadii(self):
import sys, math
minv = 0
maxv = sys.maxint * -1
for c in self.circles:
minv = min(minv, c.value)
maxv = max(maxv, c.value)
for c in self.circles:
c.r = math.pow((c.value - minv)/(maxv-minv), 0.5)*20
def project(self, globe):
# create view..
self.globe = globe
bbox = gisutils.Bounds2D()
for circle in self.circles:
x,y = globe.project(circle.lon, circle.lat)
bbox.update(gisutils.Point(x,y))
self.bbox = bbox
w = 700
self.view = gisutils.View(bbox, w, w*(bbox.height/bbox.width), 80)
# .. and place circles. you can use any other geo-libs as well, e.g. pyproj
for circle in self.circles:
x,y = self.view.project(globe.project(circle.lon, circle.lat))
circle.x = x
circle.y = y
def layout(self, steps=100):
for i in range(steps):
if i % 10 == 0:
self.toSVG()
self.layout_step()
def layout_step(self):
import math
pad = 0
for A in self.circles:
for B in self.circles:
if A != B:
radsq = (A.r+B.r)*(A.r+B.r)
d = A.sqdist(B)
if radsq + pad > d:
# move circles away from each other
v = Vector(B.x-A.x, B.y-A.y)
v.normalize()
m = (math.sqrt(radsq) - math.sqrt(d)) * 0.25
v.resize(m)
A._move(v.x*-1, v.y*-1)
B._move(v.x, v.y)
for C in self.circles:
C.move()
def toSVG(self):
from svgfig import SVG, canvas
w = self.view.width
h = self.view.height
svg = canvas(width='%dpx' % w, height='%dpx' % h, viewBox='0 0 %d %d' % (w, h), enable_background='new 0 0 %d %d' % (w, h), style='stroke-width:0.7pt; stroke-linejoin: round; stroke:#444; fill:#eee;')
g = SVG('g', id="gemeinden")
for circle in self.circles:
c = SVG('circle',cx=circle.x, cy=circle.y, r=circle.r)
c['data-key'] = circle.id
c['data-population'] = circle.value
g.append(c)
meta = SVG('metadata')
views = SVG('views')
view = SVG('view', padding="80", w=w, h=h)
proj = self.globe.toXML()
bbox = self.bbox
bbox = SVG('bbox', x=round(bbox.left,2), y=round(bbox.top,2), w=round(bbox.width,2), h=round(bbox.height,2))
views.append(view)
view.append(proj)
view.append(bbox)
meta.append(views)
svg.append(meta)
svg.append(g)
svg.save('cartogram.svg')
class Circle:
def __init__(self, lon, lat, id, value):
self.lon = float(lon)
self.lat = float(lat)
self.id = id
self.value = float(value)
self.dx = 0
self.dy = 0
def _move(self, x,y):
self.dx += x
self.dy += y
def move(self):
self.x += self.dx
self.y += self.dy
self.dx = 0
self.dy = 0
def __repr__(self):
return '<Circle lon=%f, lat=%f, id=%s, val=%f >'% (self.lon, self.lat, self.id, self.value)
def sqdist(self, circ):
dx = self.x - circ.x
dy = self.y - circ.y
return dx*dx + dy*dy
"""
been too lazy to code this myself, instead I took code from here
http://www.kokkugia.com/wiki/index.php5?title=Python_vector_class
"""
class Vector:
# Class properties
def __init__(self, x, y):
self.x = float(x)
self.y = float(y)
# represent as a string
def __repr__(self):
return 'Vector(%s, %s)' % (self.x, self.y)
'''
Class Methods / Behaviours
'''
def zero(self):
self.x = 0.0
self.y = 0.0
return self
def clone(self):
return kVec(self.x, self.y)
def normalize(self):
from math import sqrt
if self.x == 0 and self.y == 0:
return self
norm = float (1.0 / sqrt(self.x*self.x + self.y*self.y))
self.x *= norm
self.y *= norm
# self.z *= norm
return self
def invert(self):
self.x = -(self.x)
self.y = -(self.y)
return self
def resize(self, sizeFactor):
self.normalize
self.scale(sizeFactor)
return self
def minus(self, t):
self.x -= t.x
self.y -= t.y
# self.z -= t.z
return self
def plus(self, t):
self.x += t.x
self.y += t.y
# self.z += t.z
return self
def roundToInt(self):
self.x = int(x)
self.y = int(y)
return self
# Returns the squared length of this vector.
def lengthSquared(self):
return float((self.x*self.x) + (self.y*self.y))
# Returns the length of this vector.
def length(self):
from math import sqrt
return float(sqrt(self.x*self.x + self.y*self.y))
# Computes the dot product of this vector and vector v2
def dot(self, v2):
return (self.x * v2.x + self.y * v2.y)
# Linearly interpolates between vectors v1 and v2 and returns the result point = (1-alpha)*v1 + alpha*v2.
def interpolate(self, v2):
self.x = float((1 - alpha) * self.x + alpha * v2.x)
self.y = float((1 - alpha) * self.y + alpha * v2.y)
return kVec(self.x, self.y)
# Returns the angle in radians between this vector and the vector parameter;
# the return value is constrained to the range [0,PI].
def angle(self, v2):
from math import acos
vDot = self.dot(v2) / (self.length() * v2.length())
if vDot < -1.0 : vDot = -1.0
if vDot > 1.0 : vDot = 1.0
return float(acos(vDot))
# Limits this vector to a given size.
# NODEBOX USERS: name should change as 'size' and 'scale' are reserved words in Nodebox!
def limit(self, size):
if (self.length() > size):
self.normalize()
self.scale(size)
# Point Methods
# Returns the square of the distance between this tuple and tuple t1.
def distanceSquared(self, t1):
dx = self.x - t1.x
dy = self.y - t1.y
return (dx * dx + dy * dy)
# NODEBOX USERS: name should change as 'scale' is reserved word in Nodebox!
def scale(self, s):
self.x *= s
self.y *= s
return self
# NODEBOX USERS: name should change as 'translate' is reserved word in Nodebox!
def translate(self, vec):
self.plus(vec)
# NODEBOX USERS: name should change as 'translate' is reserved word in Nodebox!
def translate(self, vec, dist):
v = kVec.resize(vec, dist)
translation(v)
def distance(self, pt):
from math import sqrt
dx = self.x - pt.x
dy = self.y - pt.y
return float(sqrt(dx * dx + dy * dy))
@gka

This comment has been minimized.

Copy link
Owner Author

@gka gka commented Dec 5, 2011

actually, the trick is to find a good scale for the nodes. if you make them too large, the geo layout is too distorted. otherwise, if you make them too small, there's too much white space.. I manually adjusted the factor in line 30.

the mysterious "globe" is just a class that projects lon,lat coordinates into x,y.

@rachzhang

This comment has been minimized.

Copy link

@rachzhang rachzhang commented Sep 21, 2012

what are proj, gisutils? from django?

@almccon

This comment has been minimized.

Copy link

@almccon almccon commented Mar 4, 2014

I can't figure out proj or gisutils either, @rachzhang. @gka, are they part of Kartograph?

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