Re: [OSGeoJapan-discuss] PostGIS用の便利なファンクション
Hiroo Imaki
hiroo ＠ angeli.org
2010年 7月 27日 (火) 02:08:58 EDT
ダニエルさん、リストの皆様、
先日、ダニエルさんから指摘のあったポールさんのブログから刺激を受けてダムから河口までの距離を測るために、PostgreSQLのWITH
RECURSIVE文がどのように使えるのかをいろいろ試してみました。ポールさんのポスティングにあるSQL文だと、河川が分流したときに問題が発生するので、分流点に来たときに何かしらの基準で水が進む道を選ぶ必要があります。そこで、WITH
RECURSIVE文のUNION
ALL以下の部分で、線分の長さを使ってどちらか長いほう（または短いほう）の線分を選択するようにしてみることは出来ないかと考えたわけです。そのためにはORDER
BYとLIMITを使うことが出来るだろうと推測しました。
結果は見事に失敗で、単純にUNION ALLの文にORDER BYとLIMITは使えないようです（
http://wiki.postgresql.org/wiki/CTEReadme）。それでもしつこくいろいろ調べたり（
http://explainextended.com/2009/11/23/recursive-ctes-postgresql/）試したりしていたら、なんと、うまくいきましたので、そのSQL文をここに添付します。ただ、おはずかしながらなぜこの文でORDER
BY
とLIMITが使えたのかわからないので、ご存知の方がいたら教えてください。今書いたプロセスは、私のサイトでより詳しく報告しています。興味ある方は覗いてみてください。
http://www.geopacific.org/opensourcegis/postgis/with-recursive
WITH RECURSIVE walk_network AS (
SELECT * FROM (
SELECT * FROM network WHERE id = 1
order by length desc
limit 1) q
UNION ALL
SELECT * FROM (
SELECT n.*
FROM walk_network w JOIN network n
ON ST_DWithin(ST_EndPoint(w.segment),ST_StartPoint(n.segment),0.01)
order by n.length desc
limit 1 ) q2
)
SELECT *
FROM walk_network;
2010/7/23 Hiroo Imaki <hiroo ＠ angeli.org>
> Hi Daniel-san,
>
> Nice to talk to you over internet! Thanks for pointing out the posting. It
> was a straight forward task for me to use ESRI products but I am always
> trying to use PostGIS and other open source GIS product for my work. I love
> to try Paul's recursive method. I think it should work for my task since the
> "cost" for me is just a total length of selected segments. I even did not
> know the WITH RECURSIVE function with PostgreSQL. I definitely will learn
> about this and summarize in my site.
>
> I just recently watched "Karate-kid" not the newest one but old one with my
> kids. My five years old daughter sometimes imitates Mr. Miyagi and said
> "Daniel-san!". It just so cute to hear. Sorry for my silly comment.
>
> Thanks!
>
> Hiroo
>
> 2010/7/23 Daniel Kastl <daniel.kastl ＠ georepublic.de>
>
> Imaki-san,
>>
>> Even if I'm very glad about every pgRouting user and happy to see nice
>> tutorials in Japanese, I had to think about your water flow task when I read
>> this blog post of Paul Ramsey:
>>
>> http://blog.cleverelephant.ca/2010/07/network-walking-in-postgis.html
>>
>> Paul's example looks very much like a simple river network, and my first
>> thought was how could pgRouting eventually make use of this new feature of
>> PostgreSQL 8.4.
>>
>> Well, maybe you want to try it with your river network.
>> Not sure how to apply "costs" but maybe they are not important in your
>> case. Would be interesting to see if this is possible with only PostgreSQL
>> and PostGIS.
>>
>> Daniel
>>
>>
>>
>>
>> --
>> Georepublic UG & Georepublic Japan
>> eMail: daniel.kastl ＠ georepublic.de
>> Web: http://georepublic.de
>>
>>
>> <http://blog.cleverelephant.ca/2010/07/network-walking-in-postgis.html>
>>
>> 2010/7/22 Hiroo Imaki <hiroo ＠ angeli.org>
>>
>>> リストの皆様、
>>>
>>>
>>> 先日、pgRoutingを使ったダムから河口までの距離計算について質問した今木です。しばらくこの問題と格闘している間に、いくつか便利なファンクションを書いたので（といっても誰かの書いたものを取り込んだり、改善したりしたものですが）皆さんのお役に立てばと思いポストします。
>>>
>>> 一つ目は、複数の点（POINT)を線（LINESTRING）にスナップするファンクションです。これは、Paul
>>> Ramseyのブログからメインのコードを拝借してファンクションにしたものです。使い方は、
>>> SELECT * FROM points_snap2lines
>>> ('SELECT gid AS id, the_geom AS geom FROM point_table',
>>> 'SELECT gid AS id, the_geom AS geom FROM line_table',
>>> 0.0001);
>>>
>>> のようにして、最後に地図の単位でスナップの許容度を入力します。このファンクションを実行するとポイント、ラインのそれぞれのIDとスナップ後のジオメトリーが帰ってきます。
>>>
>>> 2つ目は、複数の点を使って複数の線をクリップするファンクションです。これは、
>>> http://postgis.refractions.net/pipermail/postgis-users/2007-September/017159.htmlの歩スティングにあったファンクションを書き換えたものです。このポスティングのままでは私のpostgresql8.3.3ではうまく動かなかったので、手直ししたのと、ST_DWithin()を使って多少、点が線からずれていてもクリップできるようにしたのと、基本的なループ上のエラーを直したものです。使い方は、
>>> SELECT * FROM split_lines2('SELECT gid AS id, the_geom AS geom FROM
>>> line_table',
>>> 'SELECT the_geom AS geom FROM point_table',
>>> 0.0001);
>>> のようにして、 ST_DWithinの距離パラメターを最後に入力します。
>>>
>>> 私の環境（WindowsXP,
>>> PostgreSQL8.3.3)ではうまく動きました。改良、バグなどありましたら連絡ください。使った感想などもお聞かせください。
>>>
>>> いまき
>>>
>>> CREATE OR REPLACE FUNCTION points_snap2lines (in pt_q text, in ln_q text,
>>> in torelance float4, out p_id int, out l_id int, out p_geom geometry)
>>> RETURNS SETOF RECORD AS
>>> $$
>>> DECLARE
>>> snap_q text;
>>> pointrec record;
>>> BEGIN
>>> EXECUTE 'CREATE TEMP TABLE line_tmp as '|| ln_q;
>>> EXECUTE 'CREATE TEMP TABLE point_tmp as '|| pt_q;
>>> snap_q :='
>>> SELECT
>>> pt_id,
>>> ln_id,
>>> ST_line_interpolate_point(
>>> ln_geom,
>>> ST_line_locate_point(ln_geom, pt_geom)
>>> ) AS the_geom
>>> FROM
>>> (
>>> SELECT DISTINCT ON (pt.id)
>>> ln.the_geom AS ln_geom,
>>> pt.the_geom AS pt_geom,
>>> ln.id AS ln_id,
>>> pt.id AS pt_id
>>> FROM
>>> point_tmp pt INNER JOIN
>>> line_tmp ln
>>> ON
>>> ST_DWithin(pt.the_geom, ln.the_geom, '||cast(torelance as text)
>>> ||')
>>> ORDER BY
>>> pt.id,ST_Distance(ln.the_geom, pt.the_geom)
>>> ) as foo' ;
>>> FOR pointrec in EXECUTE snap_q LOOP
>>> p_id := pointrec.pt_id;
>>> l_id := pointrec.ln_id;
>>> p_geom := pointrec.the_geom;
>>> RETURN NEXT;
>>> END LOOP;
>>> DROP TABLE line_tmp;
>>> DROP TABLE point_tmp;
>>> RETURN;
>>> END;
>>> $$
>>> LANGUAGE plpgsql;
>>>
>>>
>>>
>>> CREATE OR REPLACE FUNCTION split_lines2(in lineq text, in pointq text, in
>>> torelance float4, out lineid int, out the_geom geometry)
>>> RETURNS SETOF RECORD AS
>>> $$
>>> DECLARE
>>> linerec record;
>>> pointrec record;
>>> linepos float;
>>> start_ float;
>>> end_ float;
>>> loopqry text;
>>> BEGIN
>>> EXECUTE 'CREATE TEMP TABLE line_tmp as '|| lineq;
>>> EXECUTE 'CREATE TEMP TABLE point_tmp as '|| pointq;
>>> FOR linerec in EXECUTE 'SELECT * FROM line_tmp ORDER BY id' LOOP
>>> start_ := 0;
>>> end_ := 0;
>>> loopqry := '
>>> SELECT
>>> *, ST_line_locate_point('''||cast(linerec.geom as
>>> text)||''',geom) AS frac
>>> FROM point_tmp
>>> WHERE ST_DWithin(geom,'''||cast(linerec.geom as text)||''',
>>> '||torelance||')
>>> ORDER BY ST_line_locate_point('''||cast(linerec.geom as
>>> text)||''',geom)';
>>> FOR pointrec in EXECUTE loopqry LOOP
>>> end_ := pointrec.frac;
>>> lineid := linerec.id;
>>> the_geom := ST_line_substring(linerec.geom, start_, end_);
>>> start_ := end_;
>>> RETURN NEXT;
>>> END LOOP;
>>> lineid := linerec.id;
>>> the_geom:= ST_line_substring(linerec.geom, end_,1.0);
>>> RETURN NEXT;
>>> END LOOP;
>>> DROP TABLE line_tmp;
>>> DROP TABLE point_tmp;
>>> RETURN;
>>> END;
>>> $$
>>> LANGUAGE plpgsql;
>>>
>>>
>>>
>>> --
>>> Hiroo Imaki
>>> hiroo ＠ angeli.org
>>> http://www.geopacific.org
>>>
>>> _______________________________________________
>>> OSGeoJapan-discuss mailing list
>>> OSGeoJapan-discuss ＠ lists.osgeo.org
>>> http://lists.osgeo.org/mailman/listinfo/osgeojapan-discuss
>>>
>>>
>>
>> _______________________________________________
>> OSGeoJapan-discuss mailing list
>> OSGeoJapan-discuss ＠ lists.osgeo.org
>> http://lists.osgeo.org/mailman/listinfo/osgeojapan-discuss
>>
>>
>
>
> --
> Hiroo Imaki
> hiroo ＠ angeli.org
> http://www.geopacific.org
>
--
Hiroo Imaki
hiroo ＠ angeli.org
http://www.geopacific.org
-------------- next part --------------
HTML¤ÎźÉÕ¥Õ¥¡¥¤¥ë¤¬½üµî¤µ¤ì¤Þ¤·¤¿.
URL: http://lists.osgeo.org/pipermail/osgeojapan-discuss/attachments/20100726/78da5ba9/attachment-0001.html
OSGeoJapan-discuss メーリングリストの案内