PostGIS动态矢量切片（原理+实现）
矢量数据是包含空间几何字段的数据，矢量数据可视化就是将数据库中存储的矢量数据请求到前端，渲染成电子地图的过程。
为了满足用户对不同比例尺的地图的浏览需求，通常需要绘制多个比例尺级别的地图，根据用户需求加载对应的级别。但是显示器屏幕有限，当比例尺很大时，只能显示地图的局部，所以会对每个级别的地图进行分片，形成一个个小的正方形像素范围，称为瓦片，每张瓦片都有一个唯一的坐标(z, x, y)，z表示瓦片所在的层级，x表示瓦片的列号，y表示瓦片的行号。所有层级的瓦片形成了如图1所示的金字塔结构。
图1 瓦片金字塔
渲染时，前端会根据地图当前的缩放级别，确定z值，然后根据设定的显示窗口大小，确定需要加载z级的哪些瓦片，即确定多个（x, y），形成成多个瓦片坐标（z,x,y），然后根据瓦片坐标向后端请求瓦片中的数据，并渲染在窗口中。
1. 矢量数据存储
PostGIS是基于关系型数据库PostgreSQL开发的插件，用于在关系型数据库中支持空间数据的存储管理。文章MyBatisPlus+PostGIS实现Geometry数据的读写介绍了SpringBoot项目中整合矢量数据的方法。
本文中，我们以创建一张兴趣点（POI）的矢量表为例，SQL建表语句如下：
CREATE TABLE public.t_poi
(
id varchar(36) primary key, # ID
name varchar(255) not null, # 兴趣点名称
geom geometry(Point, 4326) null # 兴趣点空间位置
);
兴趣点数据存储在t_poi表里面。
2. 动态矢量切片原理
矢量瓦片是一份位于瓦片中的矢量数据的集合，只是在矢量瓦片中，这些矢量数据的空间几何字段不再是空间坐标，而是以瓦片左上角为原点的像素坐标。将原始的空间坐标系中的矢量数据映射到瓦片中，并将其空间坐标转换成瓦片中的像素坐标的过程称为矢量切片。
有两种典型的矢量切片策略：静态矢量切片和动态矢量切片。
静态矢量切片是预先将所有级别的矢量瓦片都切好，存储在文件系统中，前端请求的时候直接从文件系统中读取并返回，但是这种策略无法渲染实时更新的矢量数据，所以有了动态矢量切片策略，动态矢量切片在前端发起请求时触发的，如图2所示，根据前端传入的瓦片坐标(z, x, y)生成PostGIS的SQL语句，用PostGIS的矢量瓦片生成功能将存储的矢量数据进行坐标转换并编码成矢量瓦片。
图2 动态矢量瓦片请求示意图
接下来我们详细了解这个SQL语句，如下所示：
WITH mvtgeom AS (
SELECT id,
name,
ST_AsMVTGeom(
ST_Transform(position, 3857),
ST_TileEnvelope(z, x, y), extent => 512, buffer => 8) AS geom
FROM t_gas_station
WHERE position && ST_Transform(ST_TileEnvelope(z, x, y), 4326)
)
SELECT ST_AsMVT(mvtgeom.*) as mvt
FROM mvtgeom
SQL语句中涉及到两个坐标系4326和3857，其中4326是POI表中存储的矢量数据的空间坐标系（和建表语句对应），而3857是用于地图可视化的Web墨卡托平面投影坐标系，如图3所示。
图3 地理坐标系、Web墨卡托投影坐标系、与像素坐标系
SQL语句中，通过一系列的函数调用，将一张瓦片对应的矢量数据查询出来、将地理坐标转投影坐标、将投影坐标转像素坐标，然后将像素坐标表示的矢量数据编码便得到一张矢量瓦片，下面是每个函数的功能：
ST_Transform：该函数用于坐标转换，此处主要用于在基于球面的地理坐标系（EPSG代码为4326）和基于平面的投影坐标系(EPSG代码为3857)之间做转换。
ST_TileEnvelope(z, x, y)：该函数用于生成瓦片在Web墨卡托平面上的正方形范围。
ST_AsMVTGeom：该函数用于将Web墨卡托坐标系下的矢量数据投影到瓦片的像素坐标系中，即将矢量数据从平面坐标系换到像素坐标系中，并且只保留与落在瓦片中的部分（如图4中长方形矢量数据的阴影部分），extent参数是瓦片的像素宽高，buffer参数是瓦片向外扩展的像素数，如图4所示。向外扩展的原因是为了避免相邻瓦片拼在一起时，交界处的矢量数据在可视化效果上出现裂痕，如果两个相邻瓦片有buffer的重叠，则会消除裂痕。
ST_AsMVT：该函数将查询出来的要素编码成二进制格式的矢量瓦片，编码规则请参考MapBox定义的矢量瓦片标准。
图4 矢量瓦片的参数示意图
3. 动态矢量切片服务开发
下面我们给出基于SpringBoot和Mybatis开发的动态矢量切片的服务端代码：
import com.baomidou.mybatisplus.annotation.TableField;
import com.baomidou.mybatisplus.annotation.TableName;
import com.fasterxml.jackson.annotation.JsonIgnore;
import com.spring.accumulator.entity.handler.PointTypeHandler;
import lombok.Data;
import org.locationtech.jts.geom.Point;
/**
* 兴趣点PO
*
* @author wangrubin
*/
@Data
@TableName(value = "t_poi", autoResultMap = true)
public class PoiPO {
/**
* ID
*/
private Integer id;
/**
* POI名称
*/
private String name;
/**
* POI空间位置
*/
@JsonIgnore
@TableField(typeHandler = PointTypeHandler.class)
private Point geom;
}
import com.baomidou.dynamic.datasource.annotation.DS;
import com.baomidou.mybatisplus.core.mapper.BaseMapper;
import com.spring.accumulator.entity.PoiPO;
import com.spring.accumulator.model.vo.VectorTile;
import org.apache.ibatis.annotations.Mapper;
import org.apache.ibatis.annotations.Select;
/**
* POI表数据库访问层
*
* @author wangrubin
*/
@Mapper
public interface PoiMapper extends BaseMapper<PoiPO> {
@Select({"WITH mvtgeom AS (\n" +
" SELECT id, name, ST_AsMVTGeom(\n" +
" ST_Transform(geom, 3857),\n" +
" ST_TileEnvelope(#{z}, #{x}, #{y}), extent => 4096, buffer => 8) AS geom\n" +
" FROM t_region_poi\n" +
" WHERE geom && ST_Transform(ST_TileEnvelope(#{z}, #{x}, #{y}), 4326)\n" +
")\n" +
"SELECT ST_AsMVT(mvtgeom.*) as mvt\n" +
"FROM mvtgeom"})
VectorTile selectVectorTile(Integer z, Integer x, Integer y);
}
import com.spring.accumulator.dao.PoiMapper;
import com.spring.accumulator.model.vo.VectorTile;
import io.swagger.annotations.Api;
import io.swagger.annotations.ApiImplicitParam;
import io.swagger.annotations.ApiImplicitParams;
import io.swagger.annotations.ApiOperation;
import lombok.extern.slf4j.Slf4j;
import org.springframework.http.HttpHeaders;
import org.springframework.web.bind.annotation.GetMapping;
import org.springframework.web.bind.annotation.PathVariable;
import org.springframework.web.bind.annotation.RequestMapping;
import org.springframework.web.bind.annotation.RestController;
import javax.annotation.Resource;
import javax.servlet.http.HttpServletResponse;
import java.net.URLEncoder;
import java.nio.charset.StandardCharsets;
@RestController
@RequestMapping("/vector-tile")
public class VectorTileController {
@ApiOperation(value = "动态矢量切片请求")
@ApiImplicitParams(value = {
@ApiImplicitParam(name = "z", value = "缩放等级", required = true),
@ApiImplicitParam(name = "y", value = "瓦片行号", required = true),
@ApiImplicitParam(name = "x", value = "瓦片列号", required = true)
})
@GetMapping("/tile/{z}/{y}/{x}.pbf")
public void listPerson(@PathVariable Integer z,
@PathVariable Integer y,
@PathVariable Integer x,
HttpServletResponse response) {
try {
response.setContentType("application/x-protobuf");
response.setCharacterEncoding("utf-8");
// 这里URLEncoder.encode可以防止中文乱码
String encodedFileName = URLEncoder.encode(x.toString(), StandardCharsets.UTF_8.name()).replaceAll("\\+", "%20");
response.setHeader(HttpHeaders.CONTENT_DISPOSITION, "attachment;filename*=utf-8''" + encodedFileName + ".pbf");
VectorTile vectorTile = poiMapper.selectVectorTile(z, x, y);
response.getOutputStream().write(vectorTile.getMvt());
System.out.println(z + "-" + y + "-" + x + ":" + vectorTile.getMvt().length);
} catch (Exception e) {
// 重置response
log.error("文件下载失败" + e.getMessage());
throw new RuntimeException("下载文件失败", e);
}
}
@Resource
private PoiMapper poiMapper;
}
import lombok.Data;
import lombok.NoArgsConstructor;
@Data
@NoArgsConstructor
public class VectorTile {
byte[] mvt;
}
4. QGIS演示
在本地启动服务，端口设为8080。QGIS提供了利用矢量瓦片来渲染电子地图的功能，如图5所示，配置url为：http://localhost:8080/vector-tile/tile/{z}/{y}/{x}.pbf
图5 QGIS请求矢量瓦片渲染电子地图的配置窗口
