北京市建筑物shp数据处理：坐标校验、面积重算与人口密度估算

简介北京市城市与农村建筑物面数据采用shp格式适用于内涝治理、SWMM建模、城市规划与建筑能源分析等方向为研究者提供屋顶轮廓、建筑面积及人口分布等基础信息。压缩包共6个文件包括shp空间要素、dbf属性表、prj投影参数、shx索引及配套元数据整体约136.15MB结构清晰可直接在ArcGIS、QGIS等平台加载使用。数据字段覆盖建筑面要素、面积与人口信息可支撑下垫面分类、汇水区划分、内涝风险评估、城乡发展对比等具体工作。目前已有764人学习浏览适合需要快速获取精细建筑物轮廓数据的GIS开发者、规划师及科研人员直接投入建模分析节省数据整理与配准时间。1. 建筑物面数据与shp格式从一栋房子到一张可算的表做城市更新前期评估时经常要回答一个问题这片区域内有多少栋房子、每栋占多大面积、大概住着多少人。这类需求落到数据上最常见的就是一份北京市建筑物面数据城市农村shp格式属性表里带着建筑面要素、面积字段和人口信息。shp格式是GIS里最通用的矢量交换格式几乎任何平台都能打开但“能打开”和“能算准”是两回事——面积字段的坐标系口径、人口字段的分摊方式都会直接影响分析结论。这篇文章按我实际处理的路径来写先拆三层结构再做坐标与面积核对然后处理人口字段最后讲几个用了才知道的坑。适合刚拿到shp数据不知道怎么下手的规划从业者也适合想用建筑物数据做人口密度分析的开发同学。2. 拆开shp里的三层信息建筑面要素、面积字段与人口属性2.1 建筑面要素Polygon在shp里怎么组织建筑设计面数据里的每一个要素就是地上的一栋房子。在shp里这些房子的外轮廓以Polygon类型存储边界的每一个拐点都有明确的经度、纬度坐标。打开要素几何后看到的是一串有序坐标对首尾两点必须相同才能构成闭合的面。import geopandas as gpd gdf gpd.read_file(北京市建筑物面数据.shp, encodingutf-8) print(gdf.shape) # 行数和列数 print(gdf.geom_type.value_counts()) # 几何类型分布这一段代码做两件事第一行输出数据的行列数行数就是建筑物数量列数就是属性字段数量第二行统计几何类型正常结果应该是Polygon占绝大多数或者全部是Polygon。如果你看到很多MultiPolygon说明数据里把同一栋楼的不同楼座合并成了一个要素后续计算面积和人口关联时要注意这种多部件的特殊性。shp不是单文件而是一组同名文件共同工作。.shp存几何坐标.shx存几何索引.dbf存属性信息.prj存坐标系描述。很多人在复制数据时只拿走了.shp结果属性表全空、坐标系信息丢失这在后面会详细说。2.2 面积字段属性表里的Area为什么不能直接信拿到shp后很多人第一件事就是看属性表里的面积字段。面积字段的生成方式通常有两种一种是在ArcGIS或QGIS里用“计算几何”生成的另一种是从其他系统导入时直接带过来的。前者依赖当时图层的坐标系设定后者依赖源系统里的空间参考设置。正因为生成环境各不相同面积字段的口径差异很大。print(gdf.columns.tolist()) # 看看有哪些字段 area_col AREA # 改成实际字段名 print(gdf[area_col].describe())describe()会输出面积字段的count、mean、min、max等统计值。拿到结果先看中位数如果中位数落在几十到几百这个区间单位大概率是平方米如果中位数落在0.0001这个量级说明坐标系是经纬度几何面积是以平方度为单位算出来的这个数字没法直接用于评估。还有一种情况是面积字段有大量0值或者空值这种数据基本需要全部重算。2.3 人口信息字段语义决定数据上限这份shp里的人口信息字段是整个数据集里最需要较真的部分。建筑物本身不会“产”人口人口字段一定是从外部关联进来的。关联的源头决定了它的使用上限。常见的两种来源一种是从人口普查或社区统计按建筑物面积比例分摊到每一栋楼上的估算数另一种是房屋调查时登记的居住人数或户主登记人口属性。这两种在字段分布上有明显差异分摊得到的值通常是小数且很多栋楼的数值完全相同登记得到的值一般是整数且分布更零散。pop_col POP # 改成实际字段名 print(gdf[pop_col].describe()) print(gdf[pop_col].value_counts(dropnaFalse).head(10))第一句看人口字段整体分布第二句看重复值的头部情况。如果前10个值里有好几个人口数完全相同的记录而且数量很多基本可以判定是普查小区或居委会层级分摊到建筑上的结果。这种数据做区域汇总没有问题但精确到单栋楼就会失真。3. 用QGIS和Python核对北京建筑物shp字段、坐标系与显示设置3.1 打开shp的三种姿势与适用场景读shp数据最常见的三种工具QGIS、ArcGIS和Pythongeopandas。QGIS免费且跨平台双击拖拽就能预览适合快速目检数据范围、看属性表内容ArcGIS在企业里有大量存量模板适合要做正式制图出图的场景Python适合批量化和要写复现脚本的情况尤其是后面做字段检查、面积重算和导出。我第一次拿到这类数据时会先用QGIS拖进去看一眼再用Python做一次系统检查。纯粹用界面点来点去容易漏掉细节纯写脚本又看不到房屋和地块的空间关系两个轮着来最稳妥。打开后首先要确认能正确看到北京的轮廓如果图形出现在非洲西海岸之类的位置那基本可以判断是原始坐标系和描述信息对不上。3.2 坐标系核对CGCS2000、WGS84还是北京54北京的建筑类shp数据常见的坐标系有三种CGCS2000国家大地坐标系、WGS84GPS原始输出和北京54老项目存量数据。如果.prj里有定义Python读取后可以直接看到print(gdf.crs) # 例如 EPSG:4326 或 EPSG:4549通过EPSG编码可以快速判断坐标系类型。EPSG:4326是WGS84经纬度坐标系EPSG:4549是CGCS2000的3度带高斯投影中央经线为120度适合北京范围。规划底图通常用CGCS2000因为它是测绘成果的标准坐标框架如果从开放平台下载的数据则常见WGS84。北京54的数据只在老项目里出现拿到后最好统一转换到CGCS2000。3.3 一张自查表用5分钟确认数据能否用于后续分析我不建议拿到数据直接开算先花5分钟过一遍下面几个检查项可以省下后面排查问题的半天时间。检查项检查方法通过标准字段完整性打印列名列表有几何字段、面积字段、人口字段空间范围gdf.total_bounds输出四至经度约在[115.7, 117.4]纬度约在[39.4, 41.1]几何类型geom_type计数全部或绝大多数为Polygon面积量级面积字段describe中位数落在50至1000平方米区间人口非空人口字段isna计算空值率低于5%print(gdf.total_bounds) # 输出(minx, miny, maxx, maxy) print(gdf[pop_col].isna().mean()) # 人口字段空值率total_bounds是数据集的空间范围判断利器如果数据描述的北京却出现了经度130度、纬度20度之类的坐标说明投影信息有误或者几何本身就有问题。空值率超过20%的人口字段要先补值或重新关联不然最后的人口密度图会有大片空白区域。4. 面积重算与人口密度估算把shp变成能用的指标4.1 面积重算为什么属性表里的Area字段会翻车属性表里的面积字段不靠谱的情况我在实际项目中碰到过至少三次。表现大同小异图面上看是北京某片区面积字段却小得离谱或者单位看起来是平方米但和底图上的建筑轮廓完全对不上。原因基本是原始数据在面积计算时用了错误的坐标系或者数据经过投影转换后原有的面积字段没有同步更新。重算面积的方式不复杂先确认坐标系再投影到适合北京的等面积投影坐标系然后用几何对象的area属性重新计算# 如果读取后crs为空先要补充坐标系定义 if gdf.crs is None: gdf gdf.set_crs(EPSG:4326) # 前提是数据原本为WGS84/CGCS2000经纬度坐标 # 投影到CGCS2000 3度带中央经线120E单位变为米 gdf_proj gdf.to_crs(EPSG:4549) gdf_proj[calc_area_m2] gdf_proj.geometry.area # 用重算面积和原面积字段做个对比 gdf_proj[diff_ratio] (gdf_proj[calc_area_m2] - gdf_proj[area_col]) / gdf_proj[area_col] print(gdf_proj[diff_ratio].describe())第一段先解决crs为空的问题第二段用EPSG:4549做投影高斯-克吕格投影下的几何面积直接以平方米为单位第三段diff_ratio是重算值和原字段的差值比例。如果diff_ratio的中位数接近0说明原面积字段可信如果中位数超过0.1或出现极端值说明原面积字段失效后续统一用calc_area_m2为准。这里有一个容易被忽视的点不要用Web MercatorEPSG:3857来算面积。Web Mercator在高纬度地区面积变形严重北京的纬度高用它算出来的面积都可能偏大不少。4.2 人口字段校验同一栋楼的人口数从哪来面积重算之后下一步是校验人口字段。这一步的目的是判断人口是建筑物实测登记的还是从社区级统计往下分摊的。方法很直接看重复值、看小数位、按地块做汇总比对。# 按城市/农村分组看人口字段分布 if CATEGORY in gdf.columns: print(gdf.groupby(CATEGORY)[pop_col].describe()) # 检查人口字段是否有非整数 non_integer_mask gdf_proj[pop_col] % 1 ! 0 print(non_integer_mask.mean())如果人口字段里有大量非整数或者一组建筑的人口数完全相同基本可以判定人口是分摊出来的。分摊人口做不了单体建筑的安全评估但做街道或社区尺度的汇总分析没问题。另一套思路是用建筑面要素与人口普查单元做空间连接把普查区人口按建筑基底面积比例分摊到建筑这样得到的建筑人口在空间上更符合房屋分布但本质仍是统计推断不是实测。4.3 人口密度与导出让shp变成能进报表的指标校验完面积和人口接下来才真正构建分析指标。最常用的一个指标是人口密度但这里有一个容易混淆的点用建筑基底面积做分母得到的密度反映的是建筑覆盖区域的人口集中程度不是行政区域的人口密度。若数据里没有楼层数这个指标无法进一步折算成真实居住人口密度只能做相对比较。# 建筑人口密度单位基底面积上的人口数 gdf_proj[pop_density] gdf_proj[pop_col] / gdf_proj[calc_area_m2] # 如果城市/农村口径不同分别检查统计量 if CATEGORY in gdf.columns: print(gdf_proj.groupby(CATEGORY)[pop_density].describe()) # 导出为GeoJSON和纯表格 gdf_proj.drop(columns[pop_density]).to_file(北京市建筑物面数据_重算面积.gpkg, driverGPKG) gdf_proj[[calc_area_m2, pop_col, pop_density]].to_csv(密度指标.csv, indexFalse)到这里原始shp被加工成一张带面积、人口、密度的空间数据表可以直接用于专题制图或Excel报表。导出时选择GeoPackage而不是shp因为GeoPackage支持更长的字段名、没有.dbf的属性编码问题也避免掉shapefile本身的10字符字段名限制。5. 建筑物shp数据使用中躲不掉的5个坑从字段截断到人口口径5.1 只拷走了.shp文件属性表全空现象同事发来一份.rar解压后只有一个.shp文件用QGIS拉进去房子轮廓都在但属性表里没有面积字段也没有人口信息。原因shp的几何数据和属性数据是分文件存放的属性存在.dbf里坐标系存在.prj里只拿.shp就等于只拿到了房子的“形”丢了“身份”。解决把数据打包时确认同名的.shp、.shx、.dbf、.prj四个文件都在缺了.prj先补坐标系定义缺了.dbf基本没救只能重新从源头获取数据。5.2 字段名超过10字符直接丢内容现象用geopandas读取shp后一个字段名变成“BUILDING_”之前完整的“BUILDING_AREA_NAME”消失不见。原因DBF格式的字段名长度限制在10个字节以内中文字段名或过长英文名都会被截断截断后如果与已有字段重名还会自动加后缀。解决读进来后先打印columns发现被截断的字段先映射成可读的短名字计算完再导出导出到GeoPackage就不会再有这个烦恼。5.3 面积单位不是平方米而是度现象打印面积字段的describe结果中位数是0.0000几看建筑轮廓却明显是正常大小。原因原始数据是WGS84经纬度坐标面积计算时没投影几何面积以平方度为单位输出。解决先判断crs投影到CGCS2000 3度带或UTM 50N然后用geometry.area重算和原字段对比后修正。这里要留意转换坐标系后重算面积的操作要在投影后的图层上执行而不是在原始地理坐标图层上用面积计算工具。5.4 人口是普查单元分摊的结果现象一栋300平方米的平房和旁边一栋3000平方米的住宅楼人口字段居然都是12人。原因人口来自社区或街道级统计数据按建筑物数量平均分摊没有考虑面积和楼层差异。解决把人口字段降级为“区域汇总用”数据源在结果里标注估算口径想得到更准确的建筑人口至少要引入楼层数和建筑用途字段按“体积”分摊而不是按“数量”分摊。5.5 农村和城市的字段口径不一致现象城市建筑的人口字段有值农村区域大量为0或空值或者面积字段在城市是建筑面积、在农村是宅基地面积。原因城乡数据采集单位不同城市由房产测绘单位提供农村由乡镇或村委调查上报两套系统的字段口径没有统一。解决先查看农村与城市地块在各字段上的缺失率空值多的字段按分组单独处理统一口径时我会新建一个“数据来源级别”字段标记每个要素的口径类型防止在做区域对比时把两类数据混在一起算。6. 进阶建筑物数据的专题图与一套最小验证方案6.1 用分级着色检查城乡建筑差异数据加工完最直观的输出是做一张专题图。QGIS里右键图层选“符号化”Symbology改成“分级”Graduated字段选“pop_density”分类模式用“分位数”分级数设5档配色选从浅黄到深红。一张图上通常能看到鲜明的空间结构城区密度高、颜色深农村宅基地密度低、颜色浅。如果你想用代码验证这种分布是否合理可以分城市和农村分别看密度分位数import numpy as np if CATEGORY in gdf_proj.columns: for cat in gdf_proj[CATEGORY].unique(): subset gdf_proj[gdf_proj[CATEGORY] cat][pop_density] print(cat, np.percentile(subset, [0, 25, 50, 75, 100]))分位数输出后城市和农村的密度应该有明显分界。如果两组数据的分布完全重叠要怀疑是不是人口字段或者类别字段本身不可靠。6.2 一套最小验证方案抽样回访三件事指标做出来之后我会习惯性地抽10栋建筑做现场核对代价不大但能确认数据不是空中楼阁。抽样时优先选两类位置一类是密度最高的一档另一类是密度最低的一档。每栋楼核三件事第一建筑轮廓与现状是否一致有没有拆建后没更新的第二属性表里的人口数与实际居住规模差异是否在合理范围第三城市/农村类别标签是否准确。可以按下面的表格记录建筑编号所在区域属性表人口现场楼层数现状一致性备注001城区456一致无002农村31一致无这十分钟的抽样经常能发现“图上是空地、属性表里却有房子”之类的历史遗留问题也正好能校准密度分级的断点。能看到建筑物的实际状态后再去做密度分级我心理会踏实很多。这套路径我从第一次处理建筑shp数据用到现在希望帮到你。本文还有配套的精品资源点击获取