and the guts

Makefile

@@ -0,0 +1,38 @@
+BBOX = -122.53 37.82 -122.36 37.70
+CENTER ?= -122.44 37.758
+ZOOM ?= 13
+BBOX_FACTOR = 3.3
+
+WIDTH  ?= 3500
+HEIGHT ?= 3200
+
+# NIK2IMG ?= nik2img.py -f png -d $(WIDTH) $(HEIGHT) -b $(BBOX) --no-open
+NIK2IMG ?= nik2img.py -f png -d $(WIDTH) $(HEIGHT) -c $(CENTER) -z $(ZOOM) --bbox-factor $(BBOX_FACTOR) --no-open
+
+CHANNELS ?= trees crimes cabs
+
+all: composite.png
+
+composite.png: $(CHANNELS:%=channels/%.png)
+	python blat.py -o $@ $^
+
+all-channels: $(CHANNELS:%=channels/%.png)
+
+.SECONDARY: data/%.shp data/%.csv.vrt channels/%.png
+
+channels/%.png: data/%.shp
+	cat points.xml | perl -pi -e 's#FILENAME#$<#g' | $(NIK2IMG) - $@
+
+data/%.shp: data/%.csv.vrt
+	ogr2ogr -f "ESRI Shapefile" $@ $<
+
+data/%.csv.vrt: data/%.csv
+	cat csv.vrt | perl -pi -e 's#FILENAME#$<#g' | perl -pi -e 's#LAYER#$*#g' > $@
+
+clean:
+	rm -f *.png
+	rm -f channels/*.png
+
+data-clean:
+	rm -f data/*.csv.vrt
+	rm -f data/*.shp data/*.dbf data/*.shx data/*.prj

blat.py

@@ -0,0 +1,45 @@
+from Blit import Bitmap, Color, blends
+
+R = Color(255, 0, 0)
+G = Color(0, 255, 0)
+B = Color(0, 0, 255)
+A = Color(255, 255, 255, 0)
+WHITE = Color(255, 255, 255)
+
+def blend_raw(layers, colors, op=blends.subtract):
+    """
+    blend_raw() takes a list of Blit.Layer instances and a list of Blit.Color
+    instances, and blends them onto white with the given operation (which
+    defaults to Blit.blends.subtract).
+    """
+    out = WHITE
+    for layer, color in zip(layers, colors):
+        overlay = A.blend(color, layer)
+        out = out.blend(overlay, None, 1, op)
+    return out
+
+def blend_cmyk(red, green, blue):
+    """
+    Blend 3 layers into a single CMYK-like composite. The images are
+    interpreted as alpha masks for red, green and blue channels, which are then
+    subtractively blended onto white to simulate a CMYK print process. Remember:
+
+    white - red = cyan      (#FFF - #F00 = #0FF)
+    white - green = magenta (#FFF - #0F0 = #F0F)
+    white - blue = yellow   (#FFF - #00F = #FF0)
+
+    For best results, the layers should be grayscale with black representing
+    fully transparent and white representing fully opaque.
+    """
+    return blend_raw((red, green, blue), (R, G, B))
+
+if __name__ == "__main__":
+    import optparse
+    parser = optparse.OptionParser(usage='%prog [options] red green blue\n' + blend_cmyk.__doc__)
+    parser.add_option('--out', '-o', dest='outfile', default='blat.png')
+    options, args = parser.parse_args()
+
+    channels = map(Bitmap, args)
+    layer = blend_cmyk(*channels)
+    img = layer.image()
+    img.save(options.outfile)

csv.vrt

@@ -0,0 +1,8 @@
+<OGRVRTDataSource>
+    <OGRVRTLayer name="LAYER">
+        <SrcDataSource>FILENAME</SrcDataSource>
+        <GeometryType>wkbPoint</GeometryType>
+        <LayerSRS>WGS84</LayerSRS>
+        <GeometryField encoding="PointFromColumns" x="Longitude" y="Latitude"/>
+    </OGRVRTLayer>
+</OGRVRTDataSource>

points.xml

@@ -0,0 +1,14 @@
+<Map background-color="#000" srs="+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext  +no_defs">
+  <Style name="points">
+    <Rule>
+      <MarkersSymbolizer fill="#fff" opacity="1" width="8" height="8" stroke-width="0" allow-overlap="true" />
+    </Rule>
+  </Style>
+  <Layer name="points" srs="+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs">
+    <StyleName>points</StyleName>
+    <Datasource>
+      <Parameter name="type">shape</Parameter>
+      <Parameter name="file">FILENAME</Parameter>
+    </Datasource>
+  </Layer>
+</Map>