# -*- coding: utf-8 -*- """Convert the exported CAD building layer to GeoJSON and report what actually lands in the scene. Run with ArcGIS Desktop Python 2.7: python.exe tools/city-massing/cad-shp-to-geojson.py Only features whose envelope centre falls in the given window are written, so the GeoJSON stays a scene-sized extract of the full survey rather than a county-wide dump. """ import sys import os import io import json import arcpy def main(): if len(sys.argv) < 4: print('usage: cad-shp-to-geojson.py ') return 2 in_dir = os.path.abspath(sys.argv[1]) out_dir = os.path.abspath(sys.argv[2]) west, south, east, north = [float(v) for v in sys.argv[3].split(',')] if not os.path.isdir(out_dir): os.makedirs(out_dir) sources = [ ('Polygon-JMD.shp', 'cad-buildings.geojson', 'polygon'), ('Annotation-JMD.shp', 'cad-building-labels.geojson', 'annotation'), ('Polyline-JMD.shp', 'cad-building-outlines.geojson', 'polyline'), ] report = {'window': [west, south, east, north], 'layers': []} for filename, out_name, kind in sources: path = os.path.join(in_dir, filename) if not os.path.exists(path): print('missing %s' % path) continue fields = [f.name for f in arcpy.ListFields(path)] text_field = 'TextString' if 'TextString' in fields else None shape_type = arcpy.Describe(path).shapeType features = [] in_window = 0 total = 0 areas = [] with arcpy.da.SearchCursor(path, ['SHAPE@'] + ([text_field] if text_field else [])) as cursor: for row in cursor: total += 1 shape = row[0] if shape is None: continue extent = shape.extent cx = (extent.XMin + extent.XMax) / 2.0 cy = (extent.YMin + extent.YMax) / 2.0 if not (west <= cx <= east and south <= cy <= north): continue in_window += 1 geometry = None if shape_type == 'Polygon' and kind == 'annotation': # Annotation geometry is the label box; the point is the useful location. anchor = shape.labelPoint if hasattr(shape, 'labelPoint') else None geometry = {'type': 'Point', 'coordinates': [round(anchor.X, 7), round(anchor.Y, 7)]} if anchor else None elif shape_type == 'Polygon': rings = [] for i in range(shape.partCount): part = shape.getPart(i) ring = [[round(p.X, 7), round(p.Y, 7)] for p in part if p is not None] if len(ring) >= 4: rings.append(ring) if rings: geometry = {'type': 'Polygon', 'coordinates': [rings[0]] + [[r] for r in rings[1:]]} try: areas.append(float(shape.area)) except Exception: pass elif shape_type == 'Polyline': lines = [] for i in range(shape.partCount): part = shape.getPart(i) line = [[round(p.X, 7), round(p.Y, 7)] for p in part if p is not None] if len(line) >= 2: lines.append(line) if lines: geometry = {'type': 'MultiLineString', 'coordinates': lines} if geometry is None: continue props = {} if text_field: props['text'] = row[1] features.append({'type': 'Feature', 'properties': props, 'geometry': geometry}) out_path = os.path.join(out_dir, out_name) with io.open(out_path, 'w', encoding='utf-8') as handle: handle.write(unicode(json.dumps({'type': 'FeatureCollection', 'features': features}, ensure_ascii=False))) entry = {'source': filename, 'out': out_name, 'total': total, 'inWindow': in_window, 'written': len(features)} if areas: areas.sort() pick = lambda q: round(areas[int(q * (len(areas) - 1))], 1) entry['areaSqm'] = {'min': round(areas[0], 1), 'p25': pick(.25), 'median': pick(.5), 'p75': pick(.75), 'p95': pick(.95), 'max': round(areas[-1], 1), 'under15': sum(1 for a in areas if a < 15), 'over2000': sum(1 for a in areas if a > 2000)} report['layers'].append(entry) print('%-22s total=%6d inWindow=%6d written=%6d' % (filename, total, in_window, len(features))) if 'areaSqm' in entry: print(' area %s' % json.dumps(entry['areaSqm'])) with io.open(os.path.join(out_dir, 'cad-extract-report.json'), 'w', encoding='utf-8') as handle: handle.write(unicode(json.dumps(report, ensure_ascii=False, indent=2, sort_keys=True))) return 0 if __name__ == '__main__': sys.exit(main())