本文主要是介绍postgis计算矢量切片,希望对大家解决编程问题提供一定的参考价值,需要的开发者们随着小编来一起学习吧!
没写错,是使用postgis计算出来矢量切片。在这之前先准备一个数据:一个GIS数据表(本例中数据为一百万的点数据,坐标:4326),并在表中添加x,y字段,方便后面的数据筛选。sql中用到了
ST_AsMVT和ST_AsMVTGeom。
本文中创建矢量切片很简单,就是使用下方的一个sql,运行结果如下图。接着写一个矢量切片的http服务(参考go-vtile-example,这个例子中矢量切片压缩率更高),并且使用mapbox进行前端展示(小贴士:sql中‘points’的字符串与渲染中mapbox里的source-layer一致).代码见最下方
SELECT ST_AsMVT(tile,'points') tile FROM(
SELECT ST_AsMVTGeom(geom,ST_MakeEnvelope(100,10,125,22, 4326),4096, 0, true) AS geom FROM grid20180322 )
AS tile where tile.geom is not null
package mainimport
( _ "github.com/lib/pq""database/sql""time""log""math""errors""fmt""net/http""regexp""strconv""strings"
)func tilePathToXYZ(path string) (TileID, error) {xyzReg := regexp.MustCompile("(?P<z>[0-9]+)/(?P<x>[0-9]+)/(?P<y>[0-9]+)")matches := xyzReg.FindStringSubmatch(path)if len(matches) == 0 {return TileID{}, errors.New("Unable to parse path as tile")}x, err := strconv.ParseUint(matches[2], 10, 32)if err != nil {return TileID{}, err}y, err := strconv.ParseUint(matches[3], 10, 32)if err != nil {return TileID{}, err}z, err := strconv.ParseUint(matches[1], 10, 32)if err != nil {return TileID{}, err}return TileID{x: uint32(x), y: uint32(y), z: uint32(z)}, nil
}
type LngLat struct {lng float64lat float64
}
type TileID struct {x uint32y uint32z uint32
}func tile2lon( x int, z int)(a float64) {return float64(x) /math.Pow(2, float64(z)) * 360.0 - 180;}func tile2lat( y int, z int)(a float64) {n := math.Pi - (2.0 * math.Pi * float64(y)) / math.Pow(2, float64(z));return math.Atan(math.Sinh(n))*180/math.Pi;}func FloatToString(input_num float64) string {// to convert a float number to a stringreturn strconv.FormatFloat(input_num, 'f', 6, 64)
}func loadData(xyz TileID)(a []byte){ymax :=FloatToString(tile2lat(int(xyz.y), int(xyz.z)));ymin := FloatToString(tile2lat(int(xyz.y+1), int(xyz.z)));xmin := FloatToString(tile2lon(int(xyz.x), int(xyz.z)));xmax := FloatToString(tile2lon(int(xyz.x+1), int(xyz.z)));fmt.Println("ymax: ", ymax)fmt.Println("ymin: ",ymin)fmt.Println("xmin : ",xmin )fmt.Println("xmax : ",xmax )connStr := "dbname=xx user=xx password=xx host=localhost port=5433 sslmode=disable"db, err := sql.Open("postgres", connStr)if err != nil {panic(err)}defer db.Close()err = db.Ping()if err != nil {panic(err)}fmt.Println("Successfully connected!")var tile []bytes := []string{xmin,ymin,xmax,ymax}maxmin:=strings.Join(s, ",") s2 := []string{" where (x between", xmin,"and",xmax,") and ( y between",ymin,"and",ymax,")"}wmaxmin:=strings.Join(s2, " ") sql:="SELECT ST_AsMVT(tile,'points') tile FROM(SELECT ST_AsMVTGeom(w.geom,ST_MakeEnvelope("+maxmin+", 4326),4096, 0, true) AS geom FROM (select geom from grid20180322"+wmaxmin+") w) AS tile where tile.geom is not null"rows1:= db.QueryRow(sql)err1 := rows1.Scan(&tile)if err1 != nil {log.Fatal(err1)}fmt.Println(sql)//defer rows1.Close()return tile}
func main(){//t1 := time.Now() mux := http.NewServeMux()tileBase := "/tiles/"mux.HandleFunc(tileBase, func(w http.ResponseWriter, r *http.Request) {t2 := time.Now() log.Printf("url: %s", r.URL.Path)tilePart := r.URL.Path[len(tileBase):]fmt.Println("tilePart: ", tilePart)xyz, err := tilePathToXYZ(tilePart)fmt.Println("xyz: ", xyz)if err != nil {http.Error(w, "Invalid tile url", 400)return}tile:=loadData(xyz)// All this APi to be requests from other domains.w.Header().Set("Content-Type", "application/x-protobuf")w.Header().Set("Access-Control-Allow-Origin", "*")w.Header().Set("Access-Control-Allow-Methods", "GET, POST, OPTIONS")w.Write(tile)elapsed2 := time.Since(t2)fmt.Println("耗时: ", elapsed2)})log.Fatal(http.ListenAndServe(":8081", mux))
}
<!DOCTYPE html>
<html><head><meta charset='utf-8' /><title>Add a third party vector tile source</title><meta name='viewport' content='initial-scale=1,maximum-scale=1,user-scalable=no' /><script src='mapbox-gl.js'></script><link href='mapbox-gl.css' rel='stylesheet' /><style>body {margin: 0;padding: 0;}#map {position: absolute;top: 0;bottom: 0;width: 100%;}</style>
</head><body><div id='map'></div><script>mapboxgl.accessToken = undenfined;var tileset = 'mapbox.streets';var map = new mapboxgl.Map({container: 'map',zoom: 12,center: [109.898625809072612, 19.106708155556731],style: 'mapbox://styles/mapbox/light-v9',hash: false});map.on('load', function loaded() {map.addSource('custom-go-vector-tile-source', {type: 'vector',tiles: ['http://localhost:8081/tiles/{z}/{x}/{y}']});map.addLayer({id: 'background',type: 'background',paint: {'background-color': 'white'}});map.addLayer({"id": "custom-go-vector-tile-layer","type": "circle","source": "custom-go-vector-tile-source","source-layer": "points",paint: {'circle-radius': {stops: [[8, 0.1],[11, 0.5],[15, 3],[20, 20]]},'circle-color': '#e74c3c','circle-opacity': 1}});});</script></body></html>
这篇关于postgis计算矢量切片的文章就介绍到这儿,希望我们推荐的文章对编程师们有所帮助!