In [28]:

    

import numpy as np
import pysal as ps
import random as rdm
from pysal.contrib.viz import mapping as maps
from matplotlib.collections import LineCollection


    

In [29]:

    

shp = ps.open(ps.examples.get_path("taz.shp"))


    

In [30]:

    

dbf = ps.open(ps.examples.get_path("taz.dbf"))


    

In [31]:

    

dbf.header


    

    
Out[31]:





['AREA',
 'PERIMETER',
 'CNTY',
 'RSA',
 'AIRDB',
 'TAZ2K',
 'SQ_MILE',
 'ACRE',
 'NEWSEQ',
 'CSA',
 'CSA_NEW',
 'SQMI_TAZ',
 'TAZ_NUM',
 'CountyFIPS']

In [32]:

    

fig = figure(figsize=(9,9))
base = maps.map_poly_shp(shp)
base.set_linewidth(0.75)
base.set_facecolor('none')
base.set_edgecolor('0.8')
ax = maps.setup_ax([base])
fig.add_axes(ax)
show()


    

In [33]:

    

cnty = np.array(dbf.by_col("CNTY"))


    

In [33]:

    




    

In [34]:

    

fig = figure(figsize=(9,9))
base = maps.map_poly_shp(shp)
base.set_linewidth(0.75)
base.set_facecolor('none')
base.set_edgecolor('0.8')
counties = maps.base_choropleth_unique(maps.map_poly_shp(shp), cnty)
counties.set_linewidth(0)
counties.set_alpha(.5)
ax = maps.setup_ax([base, counties])
fig.add_axes(ax)
show()


    

In [35]:

    

cents = np.array([poly.centroid for poly in shp])
cents[0]


    

    
Out[35]:





array([  601741.78690918,  3939798.30153461])

In [36]:

    

wrook = ps.rook_from_shapefile("taz.shp")


    

In [37]:

    

w = wrook


    

In [38]:

    

cents.min(axis=0)


    

    
Out[38]:





array([  282150.8269443 ,  3615409.10372805])

In [39]:

    

def w2line_graph(w, centroids):
    
    
    
    segments = []
    for i in w.id2i:
        origin = cents[i]
        for j in w.neighbors[i]:
            dest = cents[j]
            ij = [i,j]
            ij.sort()
            segments.append([origin, dest])
    #segs = LineCollection(segments)
    
    #maps._add_axes2col(segs, [minx, miny, maxx, maxy])
    return segments


    

In [40]:

    

segs = w2line_graph(wrook, cents)


    

In [41]:

    

fig = figure(figsize=(9,9))
base = maps.map_poly_shp(shp)
base.set_linewidth(0.75)
base.set_facecolor('none')
base.set_edgecolor('0.8')
segs = LineCollection(segs)
maps._add_axes2col(segs, shp.bbox)
segs.set_linewidth(0.20)
ax = maps.setup_ax([base, segs])
fig.add_axes(ax)
show()


    

In [42]:

    

wb = ps.regime_weights(np.array(dbf.by_col("CNTY")))


    

In [43]:

    

wb.n


    

    
Out[43]:





4109

In [44]:

    

wint = ps.weights.Wsets.w_intersection(wb, wrook)


    

In [45]:

    

segs = w2line_graph(wint, cents)


    

In [45]:

    




    

In [46]:

    

fig = figure(figsize=(9,9))
base = maps.map_poly_shp(shp)
base.set_linewidth(0.75)
base.set_facecolor('none')
base.set_edgecolor('0.8')
segs = LineCollection(segs)
maps._add_axes2col(segs, shp.bbox)
segs.set_linewidth(0.20)
ax = maps.setup_ax([base, segs])
fig.add_axes(ax)
show()


    

In [47]:

    

segments = w2line_graph(wint, cents)


    

In [48]:

    

fig = figure(figsize=(9,9))
base = maps.map_poly_shp(shp)
base.set_linewidth(0.75)
base.set_facecolor('none')
base.set_edgecolor('0.8')
counties = maps.base_choropleth_unique(maps.map_poly_shp(shp), cnty)
counties.set_linewidth(0)
counties.set_alpha(.5)
segs = LineCollection(segments)
maps._add_axes2col(segs, shp.bbox)
segs.set_linewidth(0.20)
segs.set_color('0.1')
ax = maps.setup_ax([base, counties, segs])
fig.add_axes(ax)
show()


    

In [ ]: