Skip to content

Commit f66ca4e

Browse files
digoaldigoal
committed
new doc
1 parent 5d04614 commit f66ca4e

File tree

8 files changed

+192
-0
lines changed

8 files changed

+192
-0
lines changed

‎201908/20190827_01.md‎

Lines changed: 190 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,190 @@
1+
## PostGIS 3 - postgis_raster 拆离 - Separate Raster Extension
2+
3+
### 作者
4+
digoal
5+
6+
### 日期
7+
2019-08-27
8+
9+
### 标签
10+
PostgreSQL , postgis , postgis_raster , raster
11+
12+
----
13+
14+
## 背景
15+
原文
16+
17+
https://info.crunchydata.com/blog/waiting-for-postgis-3-separate-raster-extension
18+
19+
PostGIS 3开始，把raster插件单独拆出来。
20+
21+
如果要使用raster功能，需要单独编译安装
22+
23+
```
24+
CREATE EXTENSION postgis;
25+
CREATE EXTENSION postgis_raster;
26+
```
27+
28+
拆出来的好处是对于那些基本的gis用户，不需要用到raster时，不需要编译这部分功能。
29+
30+
raster拆出来，并没有减少raster的功能。
31+
32+
## 例子
33+
海平面升高30米，哪些地方会被淹掉。
34+
35+
![pic](20190827_01_pic_001.jpg)
36+
37+
创建插件，下载raster数据，并倒入数据库。
38+
39+
```
40+
# Make a new working database and enable postgis + raster
41+
createdb yvr_raster
42+
psql -c 'CREATE EXTENSION postgis' yvr_raster
43+
psql -c 'CREATE EXTENSION postgis_raster' yvr_raster
44+
45+
# BC open data from https://pub.data.gov.bc.ca/datasets/175624/
46+
# Download the DEM file, unzip and load into the database
47+
wget https://pub.data.gov.bc.ca/datasets/175624/92g/092g06_e.dem.zip
48+
unzip 092g06_e.dem.zip
49+
50+
# 倒入，并设置瓦片数量
51+
# Options: create index, add filename to table, SRID is 4269, use 56x56 chip size
52+
raster2pgsql -I -F -s 4269 -t 56x56 092g06_e.dem dem092g06e | psql yvr_raster
53+
```
54+
55+
```
56+
RELEASE: 3.0.0alpha4 GDAL_VERSION=23 (r17702)
57+
USAGE: raster2pgsql [<options>] <raster>[ <raster>[ ...]] [[<schema>.]<table>]
58+
Multiple rasters can also be specified using wildcards (*,?).
59+
60+
61+
62+
-I Create a GIST spatial index on the raster column. The ANALYZE
63+
command will automatically be issued for the created index.
64+
65+
-F Add a column with the filename of the raster.
66+
67+
-s <srid> Set the SRID field. Defaults to 0. If SRID not
68+
provided or is 0, raster's metadata will be checked to
69+
determine an appropriate SRID.
70+
71+
-t <tile size> Cut raster into tiles to be inserted one per
72+
table row. <tile size> is expressed as WIDTHxHEIGHT.
73+
<tile size> can also be "auto" to allow the loader to compute
74+
an appropriate tile size using the first raster and applied to
75+
all rasters.
76+
```
77+
78+
raster数据被分割成56x56大小的瓦片。存储在表 dem092g06e 里头，如图
79+
80+
![pic](20190827_01_pic_002.jpg)
81+
82+
如果海平面上升30米，那些瓦片会被淹掉。
83+
84+
```
85+
CREATE TABLE poly_30 AS
86+
SELECT (
87+
ST_DumpAsPolygons( -- 导出剩余空间对象（即海拔30以及以内的空间）
88+
ST_SetBandNoDataValue( -- 剔除0的数据（即海拔高于30的区域）
89+
ST_Reclass ( ST_Union(rast), -- 合并瓦片。设置0-30海平面的数据的值为1，其他为0。
90+
'0-30:1-1, 31-5000:0-0', '2BUI'),
91+
0))).*
92+
FROM dem092g06e d
93+
```
94+
95+
There are a lot of nested functions, so reading from the innermost, we:
96+
97+
- union all the chips together into one big raster
98+
- reclassify all values from 0-30 to 1, and all higher values to 0
99+
- set the "nodata" value to 0, we don't care about things that are above our threshold
100+
- create a vector polygon for each value in the raster (there only one value: "1")
101+
102+
经过以上处理，得到如下结果
103+
104+
![pic](20190827_01_pic_003.jpg)
105+
106+
有哪些建筑物被淹没了？
107+
108+
下载建筑物空间对象，倒入到buildings表
109+
110+
```
111+
# Vancouver open data
112+
# https://data.vancouver.ca/datacatalogue/buildingFootprints.htm
113+
wget ftp://webftp.vancouver.ca/OpenData/shape/building_footprints_2009_shp.zip
114+
unzip building_footprints_2009_shp.zip
115+
116+
# Options: create index, SRID is 26910, use dump format
117+
shp2pgsql -I -s 26910 -D building_footprints_2009 buildings | psql yvr_raster
118+
```
119+
120+
121+
```
122+
RELEASE: 3.0.0alpha4 (r17702)
123+
USAGE: shp2pgsql [<options>] <shapefile> [[<schema>.]<table>]
124+
OPTIONS:
125+
126+
-I Create a spatial index on the geocolumn.
127+
128+
-s [<from>:]<srid> Set the SRID field. Defaults to 0.
129+
Optionally reprojects from given SRID.
130+
131+
-D Use postgresql dump format (defaults to SQL insert statements).
132+
```
133+
134+
135+
将建筑物的srid转换为与之前的raster相同。
136+
137+
```
138+
Before we can compare the buildings to the flood zone, we need to put them into the same projection as the flood zone (SRID 4269).
139+
140+
ALTER TABLE buildings
141+
ALTER COLUMN geom
142+
TYPE geometry(MultiPolygonZ, 4269)
143+
USING ST_Transform(geom, 4269)
144+
```
145+
146+
建筑物如下
147+
148+
![pic](20190827_01_pic_004.jpg)
149+
150+
### 找到被淹建筑物的方法1，使用被淹瓦片polygon JOIN builds
151+
152+
Now we can find flooded buildings.
153+
154+
```
155+
CREATE TABLE buildings_30_poly AS
156+
SELECT b.*
157+
FROM buildings b
158+
JOIN poly_30 p
159+
ON ST_Intersects(p.geom, b.geom)
160+
WHERE ST_Intersects(p.geom, ST_Centroid(b.geom))
161+
```
162+
163+
### 找到被淹建筑物的方法2，使用直接计算被淹瓦片并 JOIN builds
164+
165+
```
166+
CREATE TABLE buildings_30_rast AS
167+
SELECT b.*
168+
FROM buildings b
169+
JOIN dem092g06e d
170+
ON ST_Intersects(b.geom, d.rast) -- 瓦片包含建筑物中心点
171+
WHERE ST_Value(d.rast, ST_Centroid(b.geom)) < 30; -- 提取瓦片内建筑物位置的海拔，当小于30时，说明这个建筑物会被淹掉。
172+
```
173+
174+
显然，方法2更加便捷，不需要落地。
175+
176+
## 参考
177+
https://info.crunchydata.com/blog/waiting-for-postgis-3-separate-raster-extension
178+
179+
http://postgis.net/docs/manual-dev/RT_ST_Value.html
180+
181+
182+
183+
<a rel="nofollow" href="http://info.flagcounter.com/h9V1" ><img src="http://s03.flagcounter.com/count/h9V1/bg_FFFFFF/txt_000000/border_CCCCCC/columns_2/maxflags_12/viewers_0/labels_0/pageviews_0/flags_0/" alt="Flag Counter" border="0" ></a>
184+
185+
186+
## [digoal's 大量PostgreSQL文章入口](https://github.com/digoal/blog/blob/master/README.md "22709685feb7cab07d30f30387f0a9ae")
187+
188+
189+
## [免费领取阿里云RDS PostgreSQL实例、ECS虚拟机](https://free.aliyun.com/ "57258f76c37864c6e6d23383d05714ea")
190+

‎201908/20190827_01_pic_001.jpg‎

42.6 KB
Loading

‎201908/20190827_01_pic_002.jpg‎

77.6 KB
Loading

‎201908/20190827_01_pic_003.jpg‎

55 KB
Loading

‎201908/20190827_01_pic_004.jpg‎

127 KB
Loading

‎201908/20190827_01_pic_005.jpg‎

116 KB
Loading

‎201908/readme.md‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -2,6 +2,7 @@
22

33
### 文章列表
44
----
5+
##### 20190827_01.md [《PostGIS 3 - postgis_raster 拆离 - Separate Raster Extension》](20190827_01.md)
56
##### 20190826_01.md [《PostgreSQL 12 + PostGIS 3 - 让空间并行计算更加智能》](20190826_01.md)
67
##### 20190825_01.md [《PostgreSQL 动态inline code - 动态存储过程代码调用 - 客户端存储过程版本管理》](20190825_01.md)
78
##### 20190821_01.md [《PostgreSQL 克隆schema - clone schema》](20190821_01.md)

‎README.md‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -54,6 +54,7 @@ digoal's|PostgreSQL|文章|归类
5454

5555
### 所有文档如下
5656
----
57+
##### 201908/20190827_01.md [《PostGIS 3 - postgis_raster 拆离 - Separate Raster Extension》](201908/20190827_01.md)
5758
##### 201908/20190826_01.md [《PostgreSQL 12 + PostGIS 3 - 让空间并行计算更加智能》](201908/20190826_01.md)
5859
##### 201908/20190825_01.md [《PostgreSQL 动态inline code - 动态存储过程代码调用 - 客户端存储过程版本管理》](201908/20190825_01.md)
5960
##### 201908/20190821_01.md [《PostgreSQL 克隆schema - clone schema》](201908/20190821_01.md)

0 commit comments

Comments
 (0)