PostGISメモ

このページでは自分がよく使うSQL文をメモしています。


SQL文をファイルから実行する

psql -U postgres -f ファイル名.sql

このコマンドを実行するためには環境変数の設定を済ませておく必要があります。


PostGISのデータベースを作成する

CREATE DATABASE sample;
\c sample
CREATE EXTENSION postgis;
\l

1行目 sampleというデータベースを作成しています。

2行目 sampleに接続します。

3行目 現在のデータベースにPostGISの拡張をインストールしています。

4行目 データベースの一覧を表示しています。


点群のテキストを読み込む

CREATE TEMP TABLE tmp
(
x double precision,
y double precision,
z double precision
);

\COPY tmp FROM 'points.txt' WITH CSV DELIMITER ' ';

CREATE TABLE pnt (
geom GEOMETRY(POINT, 4326),
elev double precision
);

INSERT INTO pnt (geom, elev) SELECT CONCAT ('SRID=4326;POINT(', x, ' ', y, ')') AS concated, z FROM tmp;

1行目 tmpという一時的に使用するデータベースを作成しています。接続が切れると自動的に削除されます。

8行目 クライアント側に存在するファイル「points.txt」をテーブル「tmp」にコピーしています。テキストはスペース区切りです。

10行目 オブジェクトタイプがポイントのpntというテーブルを作成しています。

15行目 テーブル「tmp」のフィールド「x」の値xxxxと、フィールド「y」の値yyyyを使って、文字列「SRID=4326;POINT(xxxx yyyy)」を作成し、concatedに入れています。

テーブル「pnt」のフィールド「geom」にconcatedの値を、フィールド「elev」にテーブル「tmp」のフィールド「z」の値を入れています。


ポリゴンに重なるポイントを取得する

CREATE TABLE intpnt (LIKE pnt);
INSERT INTO intpnt (geom, elev) (SELECT geom, elev  FROM pnt WHERE ST_Intersects(geom, (SELECT geom FROM mesh where コード='53381465')));

1行目 構造はpntと同じintpntというデータベースを作成しています。

2行目 テーブル「intpnt」に、テーブル「mesh」のフィールド「コード」の値が53381465のレコードのオブジェクトに重なるテーブル「pnt」のオブジェクトを追加します。


テーブル名を変更する

ALTER TABLE 〇〇〇 RENAME TO △△△;

間違えてテーブル名の先頭に数字を使ってしまったとき、〇〇〇はダブルクォーテーションで囲まないと、syntax errorになります。

PostGISのデータディレクトリを変更する


PostGISのインストール途中で「Failed to load SQL modules into the database cluster.」や「Spatial database creation failed.」とエラーが出る場合があり、データディレクトリを変更する必要がありましたので、変更手順を紹介します。


PostgreSQLのサービスを削除する

f:id:Chiakikun:20190815233634p:plain

先に設定していたデータディレクトリを削除してから、PostgreSQLインストーラーが設定したサービスを削除します。コントロールパネルから「システムとセキュリティ」→「管理ツール」→「サービス」と移動して「Postgresql-x64-10」が確認できます。

 

f:id:Chiakikun:20190815234039p:plain

cmdを管理者権限で起動します。PowerShellではダメでした。

コマンドプロンプトで次のコマンドを実行します。

sc delete postgresql-x64-10

 


データディレクトリを作成する

f:id:Chiakikun:20190815234342p:plain

データディレクトリを手動で作成します。ここはPowerShellでも構いません。initdbコマンドを使うのでPath設定はしておいてください。次のコマンドを実行します。

initdb -D C:\users\〇〇\PostgreSQL -W -U postgres -A password --no-locale -E UTF8

〇〇の部分は自身のユーザー名を入れてください。データベースのパスワードの入力を求められるのでパスワードを2回入力すると、データディレクトリが作成されます。


PostgreSQLのサービスを登録する

f:id:Chiakikun:20190815235038p:plain

サービスを登録します。cmdを管理者権限で作成して、次のコマンドを実行します。

sc create postgresql-x64-10 binPath= "\"C:\Program Files\PostgreSQL\10\bin\pg_ctl.exe\" runservice -N \"postgresql-x64-10\" -D \"C:\Users\me\PostgreSQL\"" start=auto password=************

NumbaでPythonプログラム高速化

幌延地域を対象とした10mグリッド数値標高モデルを用いた精密地形解析図の作成』を参考に、標高数値モデルから斜面方位図を作成するためのプログラムをPythonで作成したのですが、あまりにも実行速度が遅かったので、プログラムの実行速度を改善するためにNumbaというライブラリを使用し、高速化できました。Numbaについては

Pythonのnumpyで四元数してみたでも少しだけ触れたのですが、もう一度紹介します。

インストール

f:id:Chiakikun:20190825152553p:plain

管理者権限でPowershellを起動して、次のコマンドを実行します。

pip install numba

インストールはこれで終わりです。

プログラムにNumbaを適用する

①例えば、下記のコードの赤字の関数を高速化するには...

        for y in range(1, nrow - 1) :
            for x in range(1, ncol - 1) :
                if dem_check4(self.__dem, x, y, self.__nodata) == False: continue
                if x == 1:
                    dx, dy = self.__dxdy(grs80, x, y)
                dem_hoi(self.__dem, x, y, dx, dy, result)

②numbaをインポートします

# -*- coding: utf-8 -*-
import numpy as np 
import numba
....

③各関数の先頭に『@numba.jit』を記述します

@numba.jit
def dem_check4(dem, x, y, nodata):
    if dem[y][x-1]   == nodata: return False
    if dem[y][x+1]   == nodata: return False
    ....
    
@numba.jit
def dem_hoi(dem, x, y, dx, dy, result):
    Hx = (dem[y-1][x-1] + dem[y][x-1] + dem[y+1][x-1] - (dem[y-1][x+1] + dem[y][x+1] + dem[y+1][x+1])) / (3 * dx)
    Hy = (dem[y+1][x-1] + dem[y+1][x] + dem[y+1][x+1] - (dem[y-1][x-1] + dem[y-1][x] + dem[y-1][x+1])) / (3 * dy)
    ....

numbaを使ってみた結果

f:id:Chiakikun:20190825163755p:plain

f:id:Chiakikun:20190825164007p:plain

①今回は富士山周辺の数値標高モデルの斜面方位を作成しました。

numba適用前の実行時間は350秒でした。

numba適用後の実行時間は19秒でした。

ものすごく速くなったなぁと感じたのですが...

f:id:Chiakikun:20190825164702p:plain

QGISの傾斜方位を実行すると、1秒でした。

私のプログラムに何か問題があるのかもしれません...。

C# + gdal_csharpで同じプログラムを作ってみたのですが、結果は1秒切る実行速度でした。どうしても速度が必要な場合はPython以外を検討してみたほうが良いかもしれません。

最後までご覧頂き、ありがとうございました