Precomputed Gosper regions

README.md

Same as Gosper regions, but computed offline and using TopoJSON. Follow the link for an introduction.

Despite their graphical similarities, this map is 49 times bigger than the previous one (it contains 117.649 hexagonal cells), making unfeasible to render it hex by hex.

The GeoJSON hexagons are created offline by running a python script (create_hexes.py). The script is not optimized, in fact it cannot handle a 7-order Gosper curve on my machine before running out of memory.

python create_hexes.py > hexes.json

The obtained file is converted into an ESRI Shapefile with ogr2ogr, and edited with OpenJump.

ogr2ogr -f "ESRI Shapefile" hexes.shp hexes.json

By using OpenJump, Hexes with the same class are merged into a different region for each class (merge/dissolve command). The resulting shapefile is converted back to GeoJSON with ogr2ogr, then to TopoJSON, then served to the client.

ogr2ogr -f geoJSON regions.json regions.shp
topojson --cartesian --no-quantization -p class -o regions.topo.json regions.json

The rendering of TopoJSON is described in this example, while the custom projection used to make the hexes appear regular is described here. Region borders are drawn by using TopoJSON's mesh method (see the API reference and this example).

readme.regions.topo.json

create_hexes.py

# compute a Lindenmayer system given an axiom, a number of steps and rules
def fractalize(axiom, steps, rules):
    input = axiom
    for i in range(steps+1):
        output = ''
        
        for char in input:
            if char in rules:
                output += rules[char]
            else:
                output += char
                
        input = output
        
    return output
    
# convert a Lindenmayer string into an array of hexagonal coordinates
def hex_coords(fractal):
    directions = [
        {'x':+1, 'y':-1, 'z': 0},
        {'x':+1, 'y': 0, 'z':-1},
        {'x': 0, 'y':+1, 'z':-1},
        {'x':-1, 'y':+1, 'z': 0},
        {'x':-1, 'y': 0, 'z':+1},
        {'x': 0, 'y':-1, 'z':+1}
    ]
    
    # start the walk from the origin cell, facing east
    path = [{'x':0,'y':0,'z':0}]
    dir_i = 0
    
    for char in fractal:
        if char == '+':
            dir_i = (dir_i+1) % len(directions)
        elif char == '-':
            dir_i = dir_i-1
            if dir_i == -1:
                dir_i = 5
        elif char == 'F':
            dir = directions[dir_i]
            path.append({'x':path[-1]['x']+dir['x'], 'y':path[-1]['y']+dir['y'], 'z':path[-1]['z']+dir['z']})
            
    return path
    
# create a new hexagon
def new_hex(c, e):
    # conversion from hex coordinates to rect
    x = int(2*(c['x'] + c['z']/2.0))
    y = 2*c['z']
    
    return {
        'type': 'Feature',
        'geometry': {
            'type': 'Polygon',
            'coordinates': [[
                [x, y+2],
                [x+1, y+1],
                [x+1, y],
                [x, y-1],
                [x-1, y],
                [x-1, y+1],
                [x, y+2]
            ]]
        },
        'properties': {
            'class': e
        }
    }
    
# determine the order of the fractal
ORDER = 6

# create the input sequence (length equal to 7^ORDER to fill the whole fractal)
# give the element a random class (A is more likely, C is less likely)
from random import choice, randint

s = [choice(['A','A','A','A','B','B','C','D','D']) for i in xrange(7**ORDER)]

# sort the sequence by class
s.sort()

# create the Gosper curve
gosper = fractalize(
    axiom = 'A',
    steps = ORDER,
    rules = {
    'A': 'A+BF++BF-FA--FAFA-BF+',
    'B': '-FA+BFBF++BF+FA--FA-B'
})

# convert the curve into coordinates of hex cells
coords = hex_coords(fractal=gosper)

# create the GeoJSON hexes
hexes = {
    'type': 'FeatureCollection',
    'features': [new_hex(coords[i], e) for i, e in enumerate(s)]
}

# output result in JSON
import json

print json.dumps(hexes)

index.coffee

window.main = () ->
    width = 960
    height = 500
    
    svg = d3.select('body').append('svg')
        .attr('width', width)
        .attr('height', height)
        
    vis = svg.append('g')
        .attr('transform','translate(640, 120)')
        
    ### custom projection to make hexagons appear regular (y axis is also flipped) ###
    radius = 0.6
    dx = radius * 2 * Math.sin(Math.PI / 3)
    dy = radius * 1.5
    
    path_generator = d3.geo.path()
        .projection d3.geo.transform({
            point: (x,y) -> this.stream.point(x * dx / 2, -(y - (2 - (y & 1)) / 3) * dy / 2)
        })
        
    ### define a color scale (Colorbrewer's Set3) ###
    colorify = d3.scale.ordinal()
        .domain(['A','B','C','D'])
        .range(["#8dd3c7","#ffffb3","#bebada","#fb8072"])
        
    ### load topoJSON data ###
    d3.json 'readme.regions.topo.json', (error, data) ->
        ### draw the cells ###
        vis.selectAll('.region')
            .data(topojson.feature(data, data.objects.regions).features)
          .enter().append('path')
            .attr('class', 'region')
            .attr('d', path_generator)
            .attr('fill', (d) -> colorify(d.properties['class']))
            
        ### draw the boundaries ###
        vis.append('path')
            .datum(topojson.mesh(data, data.objects.regions, (a, b) -> a is b))
            .attr('d', path_generator)
            .attr('class', 'outer boundary')
            
        vis.append('path')
            .datum(topojson.mesh(data, data.objects.regions, (a, b) -> a isnt b))
            .attr('d', path_generator)
            .attr('class', 'boundary')
            

index.css

.boundary {
  stroke: #333333;
  stroke-width: 1px;
  fill: none;
}
.boundary.outer {
  stroke-width: 1.5px;
}

index.html

<!DOCTYPE html>
<html>
    <head>
        <meta charset="utf-8">
        <title>Precomputed Gosper Regions</title>
        <link type="text/css" href="index.css" rel="stylesheet"/>
        <script src="http://d3js.org/d3.v3.min.js"></script>
        <script src="http://d3js.org/topojson.v1.min.js"></script>
        <script src="index.js"></script>
    </head>
    <body onload="main()"></body>
</html>

index.js

(function() {

  window.main = function() {
    var colorify, dx, dy, height, path_generator, radius, svg, vis, width;
    width = 960;
    height = 500;
    svg = d3.select('body').append('svg').attr('width', width).attr('height', height);
    vis = svg.append('g').attr('transform', 'translate(640, 120)');
    /* custom projection to make hexagons appear regular (y axis is also flipped)
    */
    radius = 0.6;
    dx = radius * 2 * Math.sin(Math.PI / 3);
    dy = radius * 1.5;
    path_generator = d3.geo.path().projection(d3.geo.transform({
      point: function(x, y) {
        return this.stream.point(x * dx / 2, -(y - (2 - (y & 1)) / 3) * dy / 2);
      }
    }));
    /* define a color scale (Colorbrewer's Set3)
    */
    colorify = d3.scale.ordinal().domain(['A', 'B', 'C', 'D']).range(["#8dd3c7", "#ffffb3", "#bebada", "#fb8072"]);
    /* load topoJSON data
    */
    return d3.json('readme.regions.topo.json', function(error, data) {
      /* draw the cells
      */      vis.selectAll('.region').data(topojson.feature(data, data.objects.regions).features).enter().append('path').attr('class', 'region').attr('d', path_generator).attr('fill', function(d) {
        return colorify(d.properties['class']);
      });
      /* draw the boundaries
      */
      vis.append('path').datum(topojson.mesh(data, data.objects.regions, function(a, b) {
        return a === b;
      })).attr('d', path_generator).attr('class', 'outer boundary');
      return vis.append('path').datum(topojson.mesh(data, data.objects.regions, function(a, b) {
        return a !== b;
      })).attr('d', path_generator).attr('class', 'boundary');
    });
  };

}).call(this);

index.sass

.boundary
    stroke: #333
    stroke-width: 1px
    fill: none
    
    &.outer
        stroke-width: 1.5px
        

thumbnail.png