@@ -0,0 +1,1515 @@
+%22%22%22shape_msa_blockgroup.py%0A%0AOutput one shapefile per MSA containing all the blockgroups it contains%0A%22%22%22%0Aimport os%0Aimport csv%0Aimport fiona%0A%0A%0A#%0A# Import MSA to blockgroup crosswalk %0A#%0Amsa_to_bg = %7B%7D%0Awith open('data/2000/crosswalks/msa_blockgroup.csv', 'r') as source:%0A    reader = csv.reader(source, delimiter='%5Ct')%0A    reader.next()%0A    for rows in reader:%0A        msa = rows%5B0%5D%0A        bg = rows%5B1%5D%0A        if msa not in msa_to_bg:%0A            msa_to_bg%5Bmsa%5D = %5B%5D%0A        msa_to_bg%5Bmsa%5D.append(bg)%0A%0A%0A#%0A# Perform the extraction%0A#%0Afor msa in msa_to_bg:%0A    states = list(set(%5Bb%5B:2%5D for b in msa_to_bg%5Bmsa%5D%5D))%0A%0A    ## Get all blockgroups%0A    all_bg = %7B%7D%0A    for st in states:%0A        with fiona.open('data/2000/shp/state/%25s/blockgroups.shp'%25st, 'r',%0A                'ESRI Shapefile') as source:%0A            source_crs = source.crs%0A            for f in source:%0A                all_bg%5Bf%5B'properties'%5D%5B'BKGPIDFP00'%5D%5D = f%5B'geometry'%5D%0A%0A    ## blockgroups within cbsa%0A    msa_bg = %7Bbg: all_bg%5Bbg%5D for bg in msa_to_bg%5Bmsa%5D%7D%0A%0A    ## Save%0A    if not os.path.isdir('data/2000/shp/msa/%25s'%25msa):%0A        os.mkdir('data/2000/shp/msa/%25s'%25msa)%0A%0A    schema = %7B'geometry': 'Polygon',%0A              'properties': %7B'BKGPIDFP00': 'str'%7D%7D%0A    with fiona.open('data/2000/shp/msa/%25s/blockgroups.shp'%25msa, 'w', %0A            'ESRI Shapefile',%0A            crs = source_crs,%0A            schema = schema) as output:%0A        for bg in msa_bg:%0A            rec = %7B'geometry':msa_bg%5Bbg%5D, 'properties':%7B'BKGPIDFP00':bg%7D%7D%0A            output.write(rec)%0A

Add script to extract the shape of blockgroups and aggregate them per msa