Lee Hachadoorian wrote:
> I'm looking through the PostGIS reference, and I can't seem to find a
> way to take a geometry of polygons and turn it into lines. What I'm
> looking for is something like the ArcGIS Feature to Line geoprocessor,
> which will create a line shapefile where each feature is an arc
> representing the boundary between neighboring polygons with a field
> indicating the ids of the polygon on either side.
>
> Functions like ST_MakeLine require point geometries, and I don't see
> anything else that seems to be what I'm looking for. Any ideas would
> be welcome.
>
> Thanks,
> Lee Hachadoorian
> PhD Student in Geography
> Program in Earth & Environmental Sciences
> CUNY Graduate Center
This is a correction to my March 31 post. I should have have tested a little
more and not been in such a hurry I guess.
I neglected to account for the cases where a non-shared edge needs to end at
the beginning. Adding a case statement
seems to have handled that. I'd like to find a way to reference the derived
table (p4) multiple times to avoid
repeating the query (p5). The Common Table Expression WITH statement should
do this I think but it doesn't seem to be
supported in my version of PostgreSQL (8.2). Can anyone give me some ideas?
Thanks.
-Eric
/*
example_turn_poly_into_lines.sql
replace "polytable" with your poly table.
replace gid value (309) with your poly's gid.
*/
select p4.poly_id,p4.adjacent_poly_id,p4.pos_order,p4.geom
from
(
/* Get edges shared with other polys */
select t1.gid as poly_id,t2.gid as adjacent_poly_id,
ST_line_locate_point(ST_linemerge(ST_boundary(t1.geom)),ST_startpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom)))))
as pos_order,
ST_linemerge(ST_intersection(ST_linemerge(ST_boundary(t1.geom)),ST_linemerge(ST_boundary(t2.geom))))
as geom
from polytable t1, polytable t2
where t1.gid = 309
and ST_touches(t1.geom, t2.geom)
union
/*
Create edges not shared with other polys. Uses -1 for pseudo external
universe polygon.
Produces extra line(s) where there is more than one non-shared edge,
necessitating the self-join and ST_covers filter at the end.
Need help, more work on this, ...suggestions?
*/
select p1.poly_id, p1.adjacent_poly_id, p1.position as
pos_order,ST_line_substring(p3.geom,p1.position,p2.position) as geom
from
(
select t1.gid as poly_id,-1::integer as adjacent_poly_id,
ST_line_locate_point(ST_linemerge(ST_boundary(t1.geom)),ST_endpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom)))))
as position,
ST_endpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom))))
as geom
from polytable t1, polytable t2
where t1.gid = 309
and ST_touches(t1.geom, t2.geom)
and
ST_endpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom))))
not in
(
select
ST_startpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom))))
as geom
from polytable t1, polytable t2
where t1.gid = 309
and ST_touches(t1.geom, t2.geom)
)
)
as p1,
(
select t1.gid as poly_id,-1::integer as adjacent_poly_id,
case
ST_line_locate_point(ST_linemerge(ST_boundary(t1.geom)),ST_startpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom)))))
when 0 then 1
else
ST_line_locate_point(ST_linemerge(boundary(t1.geom)),ST_startpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom)))))
end
as position,
ST_startpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom))))
as geom
from polytable t1, polytable t2
where t1.gid = 309
and ST_touches(t1.geom, t2.geom)
and
ST_startpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom))))
not in
(
select
ST_endpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom))))
as geom
from polytable t1, polytable t2
where t1.gid = 309
and ST_touches(t1.geom, t2.geom)
)
)
as p2,
(
select ST_linemerge(ST_boundary(geom)) as geom
from polytable
where gid = 309
)
as p3
where p1.position < p2.position
)
as p4,
/* repeat for needed self join */
(
/* Get edges shared with other polys */
select t1.gid as poly_id,t2.gid as adjacent_poly_id,
ST_line_locate_point(ST_linemerge(ST_boundary(t1.geom)),ST_startpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom)))))
as pos_order,
ST_linemerge(ST_intersection(ST_linemerge(ST_boundary(t1.geom)),ST_linemerge(ST_boundary(t2.geom))))
as geom
from polytable t1, polytable t2
where t1.gid = 309
and ST_touches(t1.geom, t2.geom)
union
/*
Create edges not shared with other polys. Uses -1 for pseudo external
universe polygon.
Produces extra line(s) where there is more than one non-shared edge,
necessitating the self-join and ST_covers filter at the end.
Need help, more work on this, ...suggestions?
*/
select p1.poly_id, p1.adjacent_poly_id, p1.position as pos_order,
ST_line_substring(p3.geom,p1.position,p2.position) as geom
from
(
select t1.gid as poly_id,-1::integer as adjacent_poly_id,
ST_line_locate_point(ST_linemerge(ST_boundary(t1.geom)),ST_endpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom)))))
as position,
ST_endpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom))))
as geom
from polytable t1, polytable t2
where t1.gid = 309
and ST_touches(t1.geom, t2.geom)
and
ST_endpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom))))
not in
(
select
ST_startpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom))))
as geom
from polytable t1, polytable t2
where t1.gid = 309
and ST_touches(t1.geom, t2.geom)
)
)
as p1,
(
select t1.gid as poly_id,-1::integer as adjacent_poly_id,
case
ST_line_locate_point(ST_linemerge(ST_boundary(t1.geom)),ST_startpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom)))))
when 0 then 1
else
ST_line_locate_point(ST_linemerge(boundary(t1.geom)),ST_startpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom)))))
end
as position,
ST_startpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom))))
as geom
from polytable t1, polytable t2
where t1.gid = 309
and ST_touches(t1.geom, t2.geom)
and
ST_startpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom))))
not in
(
select
ST_endpoint(ST_linemerge(ST_intersection(ST_boundary(t1.geom),ST_boundary(t2.geom))))
as geom
from polytable t1, polytable t2
where t1.gid = 309
and ST_touches(t1.geom, t2.geom)
)
)
as p2,
(
select ST_linemerge(ST_boundary(geom)) as geom
from polytable
where gid = 309
)
as p3
where p1.position < p2.position
)
as p5
where ST_covers(p4.geom,p5.geom)
group by p4.poly_id,p4.adjacent_poly_id,p4.pos_order,p4.geom
having count(ST_covers(p4.geom,p5.geom)) = 1
order by p4.pos_order
---------
Eric Randall
County of Erie
140 W 6th St
Room 119
Erie, PA 16501
ph. 814-451-6063
fx. 814-451-7000
_______________________________________________
postgis-users mailing list
[email protected]
http://postgis.refractions.net/mailman/listinfo/postgis-users
