楕円内の一様乱数
などと苦労しているのですが、結局のところ、私は地図のある領域を指定して、任意の緯度・軽度情報を出す乱数を作りたかったのです。
「川の中から人が歩き始める」とか「山の中で人が消える」とか、シミュレーションと言えども、ちょっと許されない設定だと思いまして。
で、上記の記事を読まれた師匠のSさんから「こんなのがあるよ」と教え貰いました。
https://aginfo.cgk.affrc.go.jp/docs/pgisman/2.3.0/ST_GeneratePoints.html
適当なdbに接続して、例題を試してみました。
yama_db=# SELECT ST_GeneratePoints(ST_Buffer(ST_GeomFromText('LINESTRING(50 50,150 150,150 50)'), 10, 'endcap=round join=round'), 12);
st_generatepoints
--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
01040000000C0000000101000000BBA9D0DF47496240F87F69424E0E5F400101000000B7F51F68D5646240F55A898F1F17504001010000003572FBDA8A5F5F400F87342089BE5E400101000000D4A38750931D5C40ED72799B370B5E4001010000006667943B23885B4066E8F671D97C5B400101000000EBA3152848526340221E78F6C3965D400101000000113BA23198BF6240C7234B9EF33D5840010100000069499F8062745540FF88E6BF407355400101000000A581F43AC300624042D083BD5A2262400101000000636CD355C19E62404AD3293D4D904B40010100000045F7B8041DAA4E405A849AF069A34F400101000000DE1370CF02E5564050D2988109CE5940
(1 row)
geom形式で出されても全然分からんので、st_asTextでラッピングしてみました。
yama_db=# SELECT st_asText(ST_GeneratePoints(ST_Buffer(ST_GeomFromText('LINESTRING(50 50,150 150,150 50)'), 10, 'endcap=round join=round'), 12));
st_astext
----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
MULTIPOINT(147.766265850576 102.733187047914,84.0318650671786 83.712367636874,149.077046025334 68.9777848138994,107.54047530747 106.78013766944,121.059921872846 120.108716631886,137.475992983887 141.067784163524,145.074876095804 96.5277972374404,92.7965422941866 103.656943557244,66.0805226475207 56.7924152255781,72.3801321221102 71.3671567226799,145.087956158666 41.5121740980702,151.631108302923 156.218459065337)
なるほど、乱数が出力されているようです。
'LINESTRING(50 50,150 150,150 50)'), 10 → 座標 (50,50)(150,150)(150,50)で繋がれた幅10の直線上に
'endcap=round join=round'), 12)); → 12個の座標乱数を作れ
という意味のようです。
さて、実際の地図で試してみました。
kitaya_db=# SELECT st_asText(ST_GeneratePoints(ST_GeomFromText('POLYGON((35.66463989558893 139.69827111644202,35.663009879798764 139.6983247606236,35.663436999453225 139.7011571734108,35.665398233838545 139.7012966482829,35.66463989558893 139.69827111644202))'),12));
これは、以下の地図の4点で取り囲まれた地区で、任意の12点を抽出しろというSQL文になっています。
POLYGON((35.66463989558893 139.69827111644202,35.663009879798764 139.6983247606236,35.663436999453225 139.7011571734108,35.665398233838545 139.7012966482829,35.66463989558893 139.69827111644202))
ポリゴン(POLYGON)は、始点と終点(上の赤字)を同じ座標として閉じなればならないようなので、注意して下さい。
このSQL文のアウトプットは、以下のようになりました。
st_astext
-----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
MULTIPOINT(35.6638134136244 139.699085401991,35.6638440750173 139.700762425247,35.6634366319309 139.699705025931,35.6644917235626 ,35.66424835050 139.69868073913379 139.70025902483, 35.664689711471 139.700525986493,35.6635000403398 139.700601350665,35.6637472356065 139.698748086462,35.6641512918098 139.699288949827,35.6643791061995 139.701118277182,35.6636240715869 139.699272976596,35.6645803781279 139.699116246391)
エクセルで座標描いて、当ててみました。
全て領域の中に入っているようです。
さて、今度は、プログラム(go言語)でのこの座標の取り出し方です。
package main
import (
"database/sql"
"fmt"
"log"
"regexp"
// "os"
_ "github.com/lib/pq"
)
func main() {
db, err := sql.Open("postgres",
"user=postgres password=password host=192.168.0.23 port=15432 dbname=kitaya_db sslmode=disable")
if err != nil {
log.Fatal("OpenError: ", err)
}
defer db.Close()
//rows, err := db.Query("select id, age, type, departure_name, departure_number, departure_lat, departure_lng, arrival_name, arrival_number, arrival_lat, arrival_lng from user_list")
rows, err := db.Query("SELECT st_asText(ST_GeneratePoints(ST_GeomFromText('POLYGON((35.66404878 139.6931452,35.66051393 139.6943828,35.65878732 139.6973512,35.658431 139.6997472,35.66067562 139.705346,35.66404467 139.706768,35.66790807 139.7049654,35.66945399 139.702109,35.66672151 139.7018011,35.66475716 139.6987517,35.66362838 139.6955941,35.66641828 139.6934209,35.66404878 139.6931452))'),12))")
if err != nil {
log.Fatal(err)
}
defer rows.Close()
var msg string
for rows.Next() {
if err := rows.Scan(&msg); err != nil {
fmt.Println(err)
}
// まずはMULTIPOINTをバラバラに分解する
arr1 := regexp.MustCompile("[() ,]").Split(msg, -1) // '('か、')'か、' 'か、","で分割する → "[中身]" という構造でまぎらわしい
for _, s := range arr1 {
fmt.Printf("%s\n", s)
}
}
}
出力結果
> go run main3.go
MULTIPOINT
35.6601418389163
139.702654481044
35.661087233013
139.694572689447
35.6617615089132
Keyword: postgis postgres エリア 領域 範囲 指定 乱数 座標 囲む
