三島市のコインランドリー周辺500mの人口を求めてみた
人口のシェープファイルを作成する
PostGISを使ってSQLだけで面積按分:準備編（QGISを使ってCSVデータとShapeデータ結合してPostgreSQLに投入） - Qiitaqiita.com
シェープファイルからSQLに変換後、postgisのテーブル(population)へ入れる
$ shp2pgsql -s 4326 population.shp population > population.sql $ psql -d dbname -f population.sql
人口を求めるクエリを発行する
select l.id, (select name from laundry b where b.id = l.id), sum(tblt000609) from laundry l inner join population p on st_intersects( p.geom, st_buffer(l.geom::geography,500)::geometry )group by l.id
重なったところの面積を按分する
select l.id, (select name from laundry b where b.id = l.id), sum(tblt000609 *( ST_AREA( ST_Intersection(p.geom , st_buffer(l.geom::geography,500)::geometry ) ) / ST_AREA( p.geom ) ) ) from laundry l inner join population p on st_intersects( p.geom, st_buffer(l.geom::geography,500)::geometry )group by l.id
複数店舗が重なった場合、按分する
人口のシェープデータをpostgisへロードするまでの過程
人口のデータについては前の記事
shp2pgsqlコマンドを使う
シェープファイルからSQLを吐き出してくれるコマンド。
populationテーブルにインポートする場合
shp2pgsql MESH05238.shp population
吐き出したSQLは
こんな感じで人口のデータは一緒に作ってくれない。人口のデータは別立てでロードしないといけないということか。
INSERT INTO "poplation" ("mesh1_id","mesh2_id","mesh3_id","mesh4_id","obj_id","key_code",geom) VALUES ('5238','77','99','4','790694','523877994','01060000000100000001030000000100000009000000D0636666E65F61401002222222AA41402CACCCCCCC5F61401002222222AA41402CACCCCCCC5F61407896666666AA41402CACCCCCCC5F6140D82AABAAAAAA4140D0636666E65F6140D82AABAAAAAA4140741B000000606140D82AABAAAAAA4140741B0000006061407896666666AA4140741B0000006061401002222222AA4140D0636666E65F61401002222222AA4140');
人口データと結合したデータを作る
以下、とても詳しいです。QGIS上でメッシュデータと人口のデータを結合してから取り込むとのことです。
PostGISを使ってSQLだけで面積按分:準備編（QGISを使ってCSVデータとShapeデータ結合してPostgreSQLに投入） - Qiita
QGIS入門しています
レイヤーはどうも回線が遅すぎて、如何にもこうにも、まぁ、いいや。次行こう。
人口が知りたい
e-Statのダウンロードページ
静岡県は一体どのメッシュコードをダウンロードしていいのかわからない？？こちらから
検索。
レイヤに追加
上記ファイル群をダウンロードしたら。そのうちの一つ、拡張子shpをQGISのレイヤへドロップする。
レイヤの順番には意味があるので、上下で隠れてしまう場合があるので注意かも。
一応表示できたのですが、ランドマーク、星で表したところ（Postgre+postgisへinsertしたデータ)は、三島市のポイントのはずなのに、なぜか山梨に表示されてしまう。これは、チョツト問題。
測地系
国土地理院のデータは世界測地系と言うらしい。Googleさんと誤差は少ないと思われるのですが。
[GIS] 日本測地系からWGS84への変換パラメータの謎を追え - ここのことはなかったことにするかも
WGS84と緯度経度は違うようだ
Geocoding - 住所から緯度経度を検索
おそらく、緯度経度を入力してしまったから問題になったのだろう。
WGS84のSRIDは4326のようで
postgis - Choosing SRID and what is its meaning? - Geographic Information Systems Stack Exchange
postgisへinsertするときには4326でinsertすればいいようです。
QGISをインストール
Postgisの使い方
ランドマークのデータをinsert文でinsert
INSERT INTO laundry (id,name,address,geom) VALUES(1,’NARAHASHI’,’NAKATAMACHI 13-14’,ST_GeomFromText('POINT(138.5512 35.84)’, 4612));
いろいろとイメージしやすい資料
mac pgsql 使う
.bash_profileの編集
ツールへのパスを通す
export PATH=$PATH:/Applications/Postgres.app/Contents/Versions/latest/bin
db作成
CentOSにPostgreSQL/PostGISの環境構築 | JURI★GIS
created test
psql -d test
test=# create extension Postgis;
CREATE EXTENSION
test=# select postgis_full_version();
postgis_full_version
--------------------------------------------------------------------------------------------------------------------------------------------------------------------------
POSTGIS="2.2.2 r14797" GEOS="3.5.0-CAPI-1.9.0 r4084" PROJ="Rel. 4.9.2, 08 September 2015" GDAL="GDAL 1.11.5, released 2016/07/01" LIBXML="2.9.4" LIBJSON="0.12.1" RASTER
(1 row)
¥q
ファイルに用意したSQLの実行
test=# ¥i test.sql
create table
create table hoge( id int8 primary key, address varchar(200), geom geometry(point, 4326) );
