📚 课程目录第2章空间数据基础

空间数据格式与转换：Shapefile、GeoTIFF、GeoJSON、KML格式解析与互转

好，咱们直接切入正题。做空间整合分析，你手里拿到的数据格式五花八门。你可能从政府网站下载到Shapefile，从遥感平台拿到GeoTIFF，从API接口拉回一堆GeoJSON，客户又甩给你一个KML文件。如果连这些格式都玩不转，后面的分析根本没法做。

我个人习惯，拿到数据第一件事不是急着看内容，而是先搞清楚它是什么格式。为什么？因为格式决定了你能用什么工具、做什么操作，以及——最容易出bug的地方在哪。

四种格式，一招搞定

说白了，这四种格式就是不同场景下的「方言」。Shapefile是ESRI的老牌格式，GeoTIFF是栅格数据的扛把子，GeoJSON是Web端的宠儿，KML则是Google Earth的专属语言。你想想看，硬把Shapefile塞进Leaflet前端，不是不行，但没必要，转成GeoJSON效率高得多。

格式	本质	主要用途	优缺点
Shapefile	矢量数据（点线面）	ArcGIS桌面分析、传统GIS交换	字段名限制10字符、至少需3个文件
GeoTIFF	栅格数据（影像/高程）	遥感分析、DEM、地图瓦片	自带地理坐标、文件通常较大
GeoJSON	矢量数据（JSON结构）	Web地图、API传输、轻量存储	人类可读、无拓扑约束、文件易膨胀
KML	矢量+样式（XML结构）	Google Earth、移动端展示	可嵌入图标/颜色、坐标系常为WGS84

嗯，这里要注意。Shapefile其实不是「一个」文件。你拷给别人时，至少得带上.shp、.shx、.dbf这三个。我曾经碰到过同事只发了.shp给我，打开全是空图层。所以每次分享数据，我都会打包成.zip，省得丢文件。

GDAL/OGR：格式转换的瑞士军刀

做格式转换，我基本只用GDAL/OGR这套命令行工具。你问为什么不用QGIS？当然可以，但批量处理时，脚本比点击按钮靠谱一百倍。

先看最基本的转换命令。把Shapefile转成GeoJSON：

ogr2ogr -f GeoJSON output.geojson input.shp

反过来，把GeoJSON转回Shapefile：

ogr2ogr -f "ESRI Shapefile" output.shp input.geojson

GeoTIFF是栅格数据，转换逻辑不一样。用GDAL把GeoTIFF转成其它格式：

gdal_translate -of GTiff input.tif output.tif

但更常见的是，你手头有GeoTIFF，想提取里面的地理边界做成矢量——那就需要矢量化：

gdal_polygonize.py input.tif -f "GeoJSON" output.geojson

我的习惯：批量转换时，我会先写个for循环，把所有.shp遍历一遍，统一转成GeoJSON。这样前端后端都能用，省得反复沟通格式问题。

for f in *.shp; do
  ogr2ogr -f GeoJSON "${f%.shp}.geojson" "$f"
done

KML转换：看似简单，坑最多

KML的转换，我得单独拿出来说。你以为从GeoJSON转KML就是换了个后缀？天真了。

ogr2ogr -f KML output.kml input.geojson

这条命令确实能跑通。但出来的KML在Google Earth里打开，坐标系可能偏了几十米。为什么？因为KML默认使用WGS84经纬度，而GeoJSON可能是投影坐标系。

解决办法是在转换前指定源坐标系和目标坐标系：

ogr2ogr -f KML -s_srs EPSG:3857 -t_srs EPSG:4326 output.kml input.geojson

我曾经有一次帮客户做地质灾害点的可视化，他用的是西安80坐标系。我直接转了KML，结果在Google Earth里所有点位都跑到了海里。后来加上坐标系重投影才搞定。所以，坐标系这事，再怎么强调都不过分。

数据格式转换的基础知识框架

聊了这么多工具用法，我把它整理成一张结构图，方便你对照查阅。

这张图把四种格式、核心工具GDAL/OGR，以及三个关键操作串在了一起。你记住，转换不是简单的「.shp改成.geojson」。坐标系、字段名、投影信息，这些东西在转换过程中稍不注意就会丢。

Python脚本：让转换自动化

光靠命令行，手动敲起来还是累。我通常会用Python脚本封装一下，尤其是处理几十个图层时。

用Python做格式转换，核心是ogr和gdal库。安装方式很简单：

pip install gdal

一个批量转GeoJSON的脚本大概长这样：

import ogr, os, glob

input_dir = "./shapefiles/"
output_dir = "./geojson/"

for shp in glob.glob(input_dir + "*.shp"):
    src_ds = ogr.Open(shp)
    src_layer = src_ds.GetLayer()
    
    out_name = os.path.basename(shp).replace(".shp", ".geojson")
    out_path = os.path.join(output_dir, out_name)
    
    driver = ogr.GetDriverByName("GeoJSON")
    dst_ds = driver.CreateDataSource(out_path)
    dst_ds.CopyLayer(src_layer, out_name.replace(".geojson", ""))
    
    dst_ds = None
    src_ds = None
    print(f"转换完成: {out_name}")

嗯，这里要注意。如果你处理的是中文路径或中文字段名，记得加上编码设置，否则会乱码。

# 在创建数据源前设置环境变量
os.environ["SHAPE_ENCODING"] = "UTF-8"

我曾经接了个项目，数据全是用GBK编码存储的。直接用默认UTF-8读，所有村名全变成问号。后来加上SHAPE_ENCODING指定编码，才把数据救回来。遇到类似情况，你很可能也会碰到同一类坑。

避坑指南：转换中常见的那些雷

做了这么多年格式转换，我踩过的坑都够写个小册子了。这里挑几个最要命的分享给你。

避坑1：Shapefile字段名截断

Shapefile的字段名最多10个字符。如果你从GeoJSON转过来，字段名叫population_density_2023，转到Shapefile里直接被砍成populati_。
解决方案：转换时指定-lco "ENCODING=UTF-8"或者预先重命名字段。

ogr2ogr -f "ESRI Shapefile" output.shp input.geojson -lco "ENCODING=UTF-8"

避坑2：KML的坐标顺序

GDAL转换KML时，默认坐标顺序是经度纬度，也就是EPSG:4326的常规顺序。但有些老旧软件或自定义系统，读入KML时按纬度经度解析。结果就是，你的点位直接飞到地球另一端。
我曾经花了一整个下午，才发现是坐标顺序反了。检查方法：在Google Earth里打开，看看点位是不是在预期位置。

坐标系重投影：你不做，后面全白搭

格式转换里最容易被忽略的就是坐标系。我见过太多人，把WGS84的数据直接转成KML，然后在Google Earth里显示效果很好——但一做面积计算，结果明显不对。

因为WGS84是经纬度坐标系，面积计算需要投影到等面积投影。

重投影的通用做法：

# 从WGS84转到阿尔伯斯等面积投影（用于面积计算）
ogr2ogr -t_srs EPSG:5070 output_albers.shp input_wgs84.shp

# 从投影坐标系转回WGS84（用于发布到Web地图）
ogr2ogr -t_srs EPSG:4326 output_wgs84.geojson input_projected.shp

我的建议是，把数据统一到同一个坐标系再开始分析。我通常用EPSG:4326做数据存储，用EPSG:3857做Web显示，分析时再转具体需要的投影。这样既保证精度，又方便共享。

快速查看坐标系信息：拿到一个数据先看它的投影信息

ogrinfo -al input.shp -so | grep "Layer SRS"

或者用GDAL：

gdalinfo input.tif | grep -i "coord"

小结：动手试一次，比看十遍都强

格式转换这件事，说白了就是「找准工具，注意细节」。我给你的建议是，自己找几个样本数据，照着上面的命令跑一遍。跑通了，你心里就有底了。遇到报错也别慌，八成是路径问题、编码问题或者坐标系问题——按照我上面列的方法排查，基本能解决。

嗯，今天就先到这儿。你自己动手试试吧。